Dynamic response of ID bosons in a trap 
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We calculate the dynamic structure factor S(q,u>) of a one-dimensional (ID) interacting Bose gas 
confined in a harmonic trap. The effective interaction depends on the strength of the confinement 
enforcing the ID motion of atoms; interaction may be further enhanced by superimposing an optical 
lattice on the trap potential. In the compressible state, we find that the smooth variation of the gas 
density around the trap center leads to softening of the singular behavior of S(q, to) at Lieb-1 mode 
compared to the behavior predicted for homogeneous ID systems. Nevertheless, the density-averaged 
response S(q, ui) remains a non-analytic function of q and u> at Lieb-1 mode in the limit of weak trap 
confinement. The exponent of the power-law non-analyticity is modified due to the inhomogeneity 
in a universal way, and thus, bears unambiguously the information about the (homogeneous) Lieb- 
Liniger model. A strong optical lattice causes formation of Mott phases. Deep in the Mott regime, 
we predict a semi-circular peak in S(q, ui) centered at the on-site repulsion energy, u = U. Similar- 
peaks of smaller amplitudes exist at multiples of U as well. We explain the suppression of the 
dynamic response with entering into the Mott regime, observed recently by D. Clement et al., Phys. 
Rev. Lett. 102, 155301 (2009), based on an /-sum rule for the Bose-Hubbard model. 
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I. INTRODUCTION 

Experiments with ensembles of trapped cold atoms re- 
juvenated a number of directions in quantum many-body 
physics, posing new problems within that seemingly well- 
established field. One of the directions deals with the 
properties of interacting bosons confined to a strongly 
anisotropic, "one-dimensional" (ID) trap. The effective 
repulsion between bosons is enhanced by making the trap 
narrowcr[l] or by imposing a periodic potential, thus con- 
straining the kinetic energy of bosons moving along the 
trap. [2] 

The increase of the effective interparticle interaction 
affects the static and dynamic characteristics of the ID 
Bose system. Modification of the one-particle momen- 
tum distribution in the case of interaction enhanced by 
periodic potential was observed in Ref. 2. A similar ex- 
periment without a periodic potential but in a tighter ID 
trap was reported in Ref. 1. The experimental data is in 
a reasonable agreement with the predictions of the inte- 
grable Lieb-Liniger model. [3] That model allows one to 
find quantitatively the single-particle distribution func- 
tion of ID bosons at any interaction strength in the ab- 
sence of optical lattice, and establish the qualitative fea- 
tures of the distribution in the presence of the lattice, at 
least in the limit of small particle density. 

The dynamic response of ID interacting bosons can 
not be derived directly from the thermodynamic Bethe 
Ansatz solution [3] of the Lieb-Liniger model. Neverthe- 
less, a significant progress has been made in finding the 
universal singularities of the dynamic response analyti- 
cally [4-7] and determining the general structure of the 



response in a broader range of energies numerically, [8] 
utilizing the ideas of algebraic Bethe Ansatz. The most 
suitable experimental method for investigating the dy- 
namic structure factor S(q,u>) is Bragg spectroscopy. [9- 
16] By a nonlinear mixing of two laser beams, it allows 
for an independent control of the wave vector q and fre- 
quency u> of the perturbation. Recent experiments[15, 16] 
have demonstrated the use of this method in ID systems 
to study effects of interaction and periodic confinement. 

Possible complications in comparing experimental data 
with theory arise from deviations of a real atomic sys- 
tem from the ideal Lieb-Liniger model. Even in most 
favorable cases it is impossible to avoid a "soft" confin- 
ing potential applied along the ID trap, which makes the 
Bose liquid inhomogeneous. An optical lattice potential, 
applied in addition to the trap confinement to enhance 
the effective interaction between particles, brings in ad- 
ditional complications. In this paper, we develop the the- 
ory of dynamic response for an inhomogeneous ID boson 
system confined by a "soft" trap potential, both in the 
presence and in the absence of an optical lattice. In the 
absence of optical lattice, the system is a compressible in- 
homogeneous Bose liquid. The soft confinement modifies 
the dynamic response as compared to the homogeneous 
ID system. For weak interaction, the density-averaged 
response S(q, ui) of the liquid has an asymmetric peak 
as a function of frequency u at the Lieb-1 mode £+(<?), 
evaluated at the maximal gas density in the trap. The 
general shape of S(q,u>) is skewed towards low frequen- 
cies and the peak shape is a power-law singularity, with 
an exponent dependent on the interaction strength. For 
stronger interaction, the maximum shifts below e+{q). 
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In the presence of an optical lattice, the atomic cloud 
may become a mixture of incompressible Mott phases in 
equilibrium with Bosc liquid. The tight-binding limit of 
the lattice is described by the Bose-Hubbard model with 
an on-site repulsion energy U and an inter-site hopping 
matrix element J. The dynamic response of a Mott phase 
is unaffected by the trap potential and is centered at fre- 
quency u> f» U. The peak has a semicircular shape, a 
q-dependent width, and a spectral wcig ht oc (J/U). Un- 
like the Mott phase, the Bose liquid phase has dynamic 
response at small frequencies ui ~ J, and this response 
does not vanish in the limit U — ► oo. Therefore, in a 
mixed state of incompressible and compressible liquids 
the largest weight in the dynamic response comes from 
the compressible liquid and is at frequencies u> ~ J. We 
analyze in greater detail the simplest case, in which one 
Mott phase (with site occupancy 1) extends nearly over 
the entire cloud, ending with compressible-liquid caps on 
both ends. 

The paper is organized as follows. In Sec. II, we give a 
qualitative description of the main effects considered in 
the paper. In Sec. Ill, we calculate the dynamic struc- 
ture factor of the homogeneous Mott phase. In Sec. IV, 
we find dynamic responses averaged over the trap in- 
homogeneity. Attention is given to the inhomogeneous 
Bose liquid state as well as to a mixed state of coexist- 
ing Mott and Bose liquid phases. In Sec. V, we discuss 
ways to measure the dynamic structure factor, the /- 
sum rules in the presence of an optical lattice, and recent 
experiments on dynamic response of trapped atomic sys- 
tems. Appendix A deals with the finite-size effects due 
to the trap confinement. In Appendix B, we study the 
effects due to averaging over an ensemble of ID systems. 
In Appendix C, we consider dynamic response at finite 
temperatures. 



II. QUALITATIVE CONSIDERATION AND 
MAIN RESULTS 

Due to the confining potential, the density of atoms 
varies along the trap. If in addition to the potential of 
the trap an optical lattice is created, then the variation 
of the density may lead to coexistence of compressible 
and incompressible phases in the trap. In this section, 
we first identify the domain of parameters where the in- 
homogeneous cold-atom system consists of no more than 
two phases. Next, we summarize our main results for the 
dynamic structure factor averaged over the trap at zero 
temperature. 



A. Phases of cold atoms subject to the trap and 
optical lattice potentials 

The presence of an optical lattice changes the disper- 
sion relation of individual bosons. The effect of interac- 
tion becomes stronger with the narrowing of the single- 



particle bands. We consider here the limit of strong lat- 
tice potential, in which narrow bands are separated by 
wide gaps. Furthermore, we assume that only the low- 
est band is populated and thus concentrate on a one- 
dimensional Bose-Hubbard model with inter-site hopping 
matrix element J and on-site repulsion energy U, 

#BH = - J£ + b\b l+l ) + | Y, ^ { b l bl X ) 

I I 

+ 5}f(a0 (1) 

i 

V(x) = e (^y, / = 0,±1,±2,... . (2) 

Here V(x) is the confinement potential, a is the optical 
lattice constant, and fi is the chemical potential, which 
depends on the total number N of atoms in the trap. 

In the following, we view N and J as control param- 
eters (J depends exponentially on the amplitude of the 
optical lattice potential). It is convenient to introduce 
the characteristic scale, 

No = J^, (3) 
V £o 

and measure N in units of No, and J in units of U. 
The meaning of N becomes clear from considering the 
J = limit: 2Nq is the largest number of atoms fitting 
in the trap without causing double occupancy of any lat- 
tice site. In some of the existing experiments with 87 Rb 
atoms, [2, 15] the harmonic trap frequency was typically 
2tt x 60 Hz, the lattice parameter a ~ 415 nm, and inter- 
action {7/(27r) ~ 0.13 — 2 kHz which corresponds to the 
value of No in the range from 7 to 30. For investigating 
the dynamics of the compressible phase, it is desirable to 
raise the value of Nq, which may be achieved by reduc- 
ing the trap frequency. In the typical conditions of an 
experiment with multiple ID traps, [2] this would allow 
to place more atoms in each of the traps. 

Let us consider first N/Nq -C 1 and see the effect of 
decreasing J/U on the state of atomic cloud in the trap. 
At a sufficiently large J/U, the dimensionless interaction 
parameter, U / Jna, is small throughout the trap (here 
n is the position-dependent atomic density). Using the 
weak-interaction limit for the n-dependence of the chem- 
ical potential, /i = U na, and the Thomas- Fermi approx- 
imation, we find for the density no at the center of the 
trap, 

noa= {lN- ) ' (4) 

For reasons which will become clear later we refer to this 
limit as the Bogoliubov limit. The crossover from the 
weak- to strong-interaction limit occurs at the periphery 
of the cloud if J is large, and moves towards the center 
with the reduction of J/U. At N <C N , the entire cloud 
enters the strong-interaction limit before the Mott phase 
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FIG. 1: States of the trapped atomic cloud in the presence 
of optical lattice. Compressible (C) and incompressible (I) 
phases coexist in distinct states, such as CIC, etc. The solid 
line serves as boundary between these states. Its first two 
segments (denoted by CD and CD) bound the state CIC, which 
consists of an incompressible region in the middle of the trap 
and two compressible ones at the sides, see Fig. 2. The in- 
compressible region nucleates at line © upon crossing over 
from state C by decreasing J/U or increasing N/Nq. Simi- 
larly, traversing line © upwards leads to a compressible region 
form in the trap center; hence the combination CICIC. The 
small- J asymptotes of lines © and © are given in Eqs. (7) and 
(57), respectively. The crossover line (dashed) at N/No < 1 
is defined by Eq. (6) and separates the regimes of weak and 
strong interaction in the compressible state. At N/No > 1, 
a mixed state with compressible and incompressible phases 
forms to the left of the vertical line (dotted). The compress- 
ible phases represent a substantial fraction of the cloud in 
the shaded region; its left boundary is given by Eq. (6). To 
the left of that line, the cloud is mostly distributed between 
incompressible states. We are interested in the evolution of 
the dynamic structure factor upon decreasing J/U at small 
N/Nq, dash-dotted line. 



nucleates. The strong-interaction limit, U/Jna 1, is 
described by the so-called Tonks-Girardeau gas, for which 
the analog of Eq. (4) reads 



2 i/ 2 ^jvy/a 

iVo 



1/4 



n a <C 1. 



(5) 



The crossover between the Bogoliubov and Tonks- 
Girardeau limits occurs at U / Jn^a ~ 1. Using the es- 
timate for n (either Eq. (4) or Eq. (5)), we find for the 
crossover line 



TV 



3/2 



(6) 



(Here we dispensed with the numerical factors present 
in Eqs. (4) and (5), as Eq. (6) defines a crossover rather 
than a transition.) The crossover line is shown in Fig. 1 
(see dashed line at N/No < 1). Further decrease of J 
(moving along the dash-dotted line in Fig. 1) eventually 
allows one to reach the boundary between the Tonks- 
Girardeau gas and the Mott phase. The nucleation of 




FIG. 2: Density profile and energy diagram of an atomic cloud 
in the Mott-insulator regime confined in a trap. Upper panel: 
Energy diagram showing two lowest Mott bands curved due 
to the trap potential. The hashed area shows the phase space 
occupied by atoms below the level of the chemical potential 
p. To leading order of perturbation theory the bands are 
separated (center to center) by energy U and have widths of 
4J and 8J for the lower and upper bands, respectively. The 
state of the liquid depicted here corresponds to the state CIC 
of Fig. 1. Lower panel: Density profile n(x) of the atomic 
cloud in the trap. The incompressible region with na = 1 
occupies the center of the trap for 2J < y, < U — 4J at 
J/U <C 1 and it borders with two compressible regions on the 
sides. The compressible regions have n(x) given in Eq. (59). 



the Mott phase occurs in the center of the trap, once the 
condition N = (8/tt)^/ J/cq is reached; in dimcnsionlcss 
variables, the latter condition reads 



N 



1/2 



(7) 



The corresponding boundary is shown in in Fig. 1 by the 
solid line (segment ©). Under the condition N/Nq <C 1, 
the Mott phase occurs [17] at J/U <C 1. The advantage 
of a small number of particles in the trap, N <C Nq, is 
that there is at most one domain of the trap occupied by 
the Mott phase, see Fig. 2. 

For a large number of atoms, N ^> Nq, the composition 
of the atomic cloud in the trap may become more com- 
plicated. Upon the reduction of the bandwidth, the first 
incompressible phase appears [19, 20] at J/U < 0.28, see 
dotted line in Fig. 1. Further reduction of the band width 
brings an increasing number of Mott phases correspond- 
ing to various integer values j of the filling, n(l)a = j. 
Incompressible phases nucleate around multiple points in 
the trap, where the condition J/U ~ 1/j is satisfied [21]. 
Once the band width is reduced to satisfy the condition 
J/U ~ (n a) _1 , a major part of the cloud belongs to 
the incompressible phases. The value of n^a is estimated 
by Eq. (4) both in the limits of weak (J/U ^> 1) and 
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strong (J/U <C (noa)" 1 ) interaction, where the equation 
of state for incompressible phases is applicable. Using 
Eq. (4), we may re-write the condition J/U = (n-oa) -1 
separating "mostly compressible" state from the "mostly 
incompressible" one in the form of Eq. (6). 

We summarize the important for us features of the 
phase diagram for the trapped cloud subject to the op- 
tical lattice potential in Fig. 1. Reducing the parameter 
J/U at N/Nq <C 1 allows one to investigate the gradual 
increase of correlations due to interaction in the liquid 
phase as well as the formation of a single Mott phase do- 
main in contact with Bose liquid (see Fig. 2). In the next 
section, we review the manifestation of such an evolution 
of the state of the cloud in the dynamic structure factor. 



B. Average dynamic structure factor in the 
absence of optical lattice 

In a homogeneous ID system with contact interaction 



g5(xi 



between bosons (Lieb-Liniger model [3]), the 



dynamic response depends on the dimcnsionless interac- 
tion constant 7 = mg/n. In this section, we analyze the 
density-averaged dynamic structure factor for a gas in a 
trap. 

The dynamic structure factor at zero temperature is 
defined as S(q,tu) = lm(p(x,t)p(0,0)) q ^, where p(x,t) 
is the density fluctuation operator. In the Lieb-Liniger 
model [3] the dynamic structure factor has a power-law 
divergence [4, 5] at the so-called [3] Licb-1 mode, to = 
£+(q), and vanishes at the threshold determined by the 
spectrum of the Lieb-2 mode: S(q,uj) = at u> < S-(q). 

The full form of S(q,uj) can be found analytically in 
the limit of infinitely-strong contact repulsion between 
the bosons (7 3> 1, limit of the Tonks- Girardeau gas): 



m 
<1 

nnq 
m 



9(fj-e-(q))0(e+(q)-u), 
2 



± 



Q 
2m 



(8) 



If 7 <C 1, then the dispersion of Lieb-1 mode approaches 
the Bogoliubov spectrum, 



with the collective velocity v = y^gn/m. 
expect, S(q,ui) approaches a <5-function, 



(9) 



As one may 



S(q,w) 



nnq 



me- 



-(?) 



-(<?)) 



(10) 



in this limit. [4] The Lieb-2 mode in the 7 <C 1 limit cor- 
responds to the grey solitons spectrum. [22] The spectral 
weight near e_ (q) is, however, extremely small [23] at 
any q. 




FIG. 3: Dynamic response S(q, to) as a function of u for an 
atomic cloud in the Tonks- Girardeau limit, see Sec. IV C. 
The non-analyticities are at w = e±(q), with e±(q) given in 
Eq. (8). The plot is made using Eq. (13) for the choice of 
parameters q — 1.2-jmo. In the presence of optical lattice, one 
replaces m —* m* . 



At arbitrary 7, the general form of S(q,uj) is [4, 5] 



S(q,u) 



5e 



LU — Ej 



[0(e+ + vi6(u - £+)] 



(11) 

with the exponent p\ and asymmetry factor v± depending 
on momentum q and density n. 

All the above forms of S(q,u>) were derived for a 
spatially-homogeneous Bose liquid. How one may com- 
pare the theory for a homogeneous system with measure- 
ments performed in the presence of a trap? We consider 
here the most favorable case of a smooth longitudinal po- 
tential confining the ID bosons, whose confinement en- 
ergy o^o is negligibly small (not resolved on measurement 
scale). Due to the interaction between bosons, the system 
length 2L is usually much larger than the quantum length 
scale A of the confinement, allowing for independent (un- 
affected by interference) probing of different parts of the 
system. Provided the system is probed at a sufficiently 
large momentum, g > 1/A, one may regard portions of 
length A ~ 1/ -y/moJo as homogeneous and responding in- 
dependent of each other. Then the dynamic structure 
factor can be approximated by the density-averaged one, 



S(q,u>) 



1 

2L 



dxS (q, oj; n(x)) . 



(12) 



The confinement affects the dynamic response through 
the variation of the density profile n(x). The averaged 
structure factor S(q, u>) may differ considerably from 
S(q,u)) of a homogeneous system. We note that the 
accuracy of the density-averaged approximation applied 
to a smooth trapping potential was checked in the limit 
7 — > 00 and yielded excellent results [24] (error less than 
5%). Constraints on the applicability of Eq. (12) may 
arise in some cases due to quantum corrections to the dy- 
namic response originating from momentum uncertainty 
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FIG. 4: Dynamic response S(q, u>) as a function of u> for an 
atomic cloud in the Bogoliubov limit, see Sec. IV B. A square- 
root singularity occurs at u = £+(?), where £+(?) is given by 
the Bogoliubov spectrum in Eq. (9). The plot is made using 



Eq. (14) for the choice of parameters q /2m 



In the 



presence of optical lattice, one replaces m — > m* and g —* Ua. 



on the scale A -1 . We defer the analysis of this "finite- 
size" effect to Appendix A. 

For dcfinitcncss, we assume here a parabolic confine- 
ment, see Eq. (2). However, our main conclusions require 
only the existence of a smooth maximum in the density 
profile, while the specific form of the n(x) dependence is 
not important. For the parabolic confinement, averaging 
of Eq. (8) over the proper density profile yields [24] 



S(q,uj) 



7m <r 



- 9 (£_(<?) - w) ^{e+(q)+uj) (e-(q)-io) 
x6(e+(q)-u>), 7 »1, (13) 

where for brevity we included only interval of wavevectors 
q < 27rno (see Sec. IV C for larger q). The w-dependcncc 
of S{q 1 Lo) in Eq. (13) for a fixed value of q is shown in 
Fig. 3. 

Similarly, averaging of Eq. (10) yields 

B( \ nn ° ^ ~ (l 2 / 2m ) 2 Q/ , \ \ 

S{q,uj) = 3 -p(e+(g) - u) 

x0(uj-q 2 /2m), 7 <1. (14) 

The dependence of S(q,u>) in Eq. (14) on u) is shown in 
Fig. 4. The parameters n and v in the dispersion rela- 
tion £+(<?), entering Eqs. (13) and (14), correspond to 
the maximal density, n = no. While the full form of the 
averaged structure factor substantially uses Eqs. (8) and 
(10) and assumes the parabolic confinement, the behav- 
ior of S(q,u>) near the upper edge, uj — > £+(<?), is uni- 
versal and relies only on the form of the corresponding 
singularity in S(q,u>) before the averaging, and on the 
existence of quadratic maximum in the function n(x). 
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FIG. 5: Sketch showing the singular behavior of S(q,ui) in 
Eq. (15) at the Lieb-1 mode. The choice of parameters corre- 
sponds to the regime of "large q and arbitrary 7" considered 
in Sec. IV D 1. The dynamic response diverges at the Lieb-1 
mode £+(<?) (solid line) and is strongly suppressed away from 
£+{q) at distances exceeding vq (dotted line). The solid line 
is plotted using Eq. (88) with the choice ^1 = 0.8; the dotted 
line is hand-drawn and illustrates the qualitative behavior of 
S(q,u}) away from the Lieb-1 mode. 



Tuning of interaction from weak to strong may be per- 
formed, for example, by decreasing the number of atoms 
in the trap, which reduces no . The response should then 
cross over from the one with divergence at the upper edge 
£+(q), see Eq. (14) and Fig. (4), to the one with a non- 
analytic point at uj = s+(g) and maximum below £+(<?), 
see Eq. (13) and Fig. (3). 

At arbitrary interaction, £'((?, u>) demonstrates a non- 
analytical behavior in the vicinity of uj = £+(q) of the 
general type 



S(q,i 



Sc 



uj - £+{q) 
+B6(u;-£ + (q))} + C, 



[A9(e + (q) - uj) 



(15) 



with a g-dependent exponent /ii evaluated at maximum 
density n = no; constants A, B, and C are given in Sec- 
tion IV D. We illustrate the singular behavior of S(q,ui) 
at fix > 1/2 and v\ = 1 in Fig. 5. The divergence in 
S(q,ui) disappears at < 1/2, but the point u — £+{q) 
remains non-analytic. 



C. The dynamic structure factor in the presence of 
optical latice, N/N a < 1 

In the above consideration, we assumed no optical 
lattice, so that the dispersion relation of ID bosons is 
quadratic. The presence of an optical lattice changes 
the dispersion relation of individual bosons; sufficiently 
strong lattice potential results in narrow bands of al- 
lowed energies. Equations (13)-(15) remain true even in 
that case, provided N < and J/U > (N/2N ) 2 , see 
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Eqs. (6) and (7), and assuming go < 1. In the narrow- 
band limit, the free-particle mass m in Eqs. (13)-(15) 
must be replaced by l/2Ja 2 , where 4J is the bandwidth 
of the lowest allowed energy band. The general form of 
the singular contribution, Eq. (15), is applicable also in 
the crossover regime specified by Eq. (6). 

The condition J/U 3> (N/2Nq) 2 ensures small occupa- 
tion of the lattice sites, n^a <C 1, which allows one the use 
of parabolic dispersion relation for bosons. A reduction of 
the bandwidth eventually leads to occupation n a > 1/2 
in the center of the trap, see Eq. (5). Under this condi- 
tion, the inflection points in the full dispersion relation 
for a narrow-band spectrum become important, leading 
to specific singularities in the dynamic structure factor. 
These singularities, associated with symmetric configura- 
tions of excited particles and holes, where studied for the 
homogeneous case in the context of spin chains, [25-27] 
and we will not dwell on them here. 

Further decrease of the bandwidth results in the for- 
mation of an incompressible Mott phase in the center of 
the trap at n a a = 1. The Mott phase turns, on both 
sides, into compressible regions. In these compressible 
caps, the concentration n drops from n = l/a to zero 
continuously, see Fig. 2. The dynamic response is, there- 
fore, equal to the sum of contributions of the Mott phase 
(incompressible) and two caps (compressible), 



S(q,uj) = S inc (q,oj) + 2S com (q,oj). 



(16) 



Here, the response is averaged over the trap with an 
appropriate density profile, see Fig. 2, which accounts 
for the presence of the compressible and incompressible 
parts. At J/U < (N/2N Q ) 2 , the Mott phase extends 
almost over the entire cloud, with the caps constitut- 
ing a small part of the system. To the leading order in 
J/U -C 1, the Mott phase contribution is given by 



Si nc (q,(j) 



64 J 2 sin 2 (<ja/2) 
aU 2 Wi(q) 



1 - 



uj-U 
Wx{q) 



(17) 



where W\(q) = 2J yf\ + 8 cos 2 (qa/2) and we assumed 
one atom per site (N -C No), see Section IIIB for de- 
tails. Thus, at the leading order of J/U, the Mott phase 
produces a semicircular peak of width 2W±(q), centered 
at u) = U. The contribution of the liquid caps at frequen- 
cies u> w U can be safely neglected since it is small by 
a factor AN/N < 1, where AN = 4J/e N is the num- 
ber of atoms in the caps. At much smaller frequencies, 
u> ~ J, the dynamic response is solely due to the caps. 
The contribution of one cap to Eq. (16) is 

= , 2uj 9(oj)d(4J S m(qa/2)~uj) 

Jcoml?,") — — — , [iS) 



ae ° N2 y/[U&m(qa/2)] 2 - ^ 



see Section IV E. Here, the divergent contribution at the 
threshold frequency comes from parts of the trap with oc- 
cupation na m 1/2 and has the origin in the particle-hole 
symmetry [25-27] mentioned above. The total dynamic 
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FIG. 6: Dynamic response S = Sine + 25 C omp as a function of 
ui for an atomic cloud in the Mott regime, featuring a square 
root singularity at frequency cj* = 4Jsin(ga/2) and a semi- 
circular peak at ui = U. The plot is made using Eqs. (16), 
(17), and (18) with the parameters values: q — n/a, U = 10J, 
and e iV 2 = 10J. 



response in the Mott regime is shown in Fig. 6 and is 
plotted using Eqs. (16), (17), and (18). 

The dynamic response S(q, oj) in the Mott regime 
(J <C U) is peculiar, featuring a hierarchy of weights 
at different frequencies. The smallest frequency scale is 
uj ~ J, where the dynamic response is independent of 
interaction at large U . The other frequencies are integer 
multiples of U , uj ~ jU, where the dynamic response van- 
ishes with increasing U as oc (J/U) 2 '' , for j = 1,2,..., 
as long as the applicability of the Bose-Hubbard model 
holds. Ignoring for simplicity the fine structure of the 
peaks (each peak width is proportional to J), we repre- 
sent the dynamic structure factor as 



S{q,uj) 



J2 sin 1 ^ 
j=o 



(qa/2) 



S(uj~n 3 (q)). 



where i>(0) = 1 and v(j > 0) 



(19) 

2, fl = 4Jsin(ga/2) and 



j>0 



jU, and ctj are numerical coefficients. 



A possible way to tune between different states of a 
ID Bose system is provided by control of the optical lat- 
tice potential. Typically [2] the crossover between the 
limits of weak and strong interaction occurs for a suffi- 
ciently strong optical lattice potential, allowing one the 
use of the tight-binding approximation for the single- 
particle spectrum and of the ID Bose-Hubbard model 
for the description of interaction. The increase of the op- 
tical lattice potential leads chiefly to the decrease of the 
single-particle band width 4 J, while the on-site repulsion 
U changes little. Narrowing of the band width 4J leads 
to the evolution of the dynamic structure factor. This 
evolution allows for a clear interpretation at N <C No, 
see Sec. II A. We will follow the changes in the dynamic 
response accompanying the reduction of the J/U ratio at 
small N, e.g., along the dash-dotted line in Fig. 1. 

The evolution between the regimes of weak interaction, 
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strong interaction in the compressible phase, and the 
regime of coexistence of compressible and Mott phases 
is best illustrated by the integral characteristics of the 
dynamic response. For the weakly-interacting gas at 
qa <C 1, we find 

[ du- f . 3^ / N\ 1/3 (J\ 1/2 
J^ S ^ = I^{-N ) [u) * (20) 

while after the crossover to the Tonks-Girardeau limit 
one obtains 

(-) 

The above two estimates match each other at the 
crossover line in Eq. (6). Remarkaby, in the Tonks- 
Girardcu limit the integral intensity of the dynamic re- 
sponse is independent of the parameters of the system. 

Upon further bandwidth reduction, the Mott phase is 
formed, and at J/U < (nN/8N ) 2 , it occupies the most 
part of the trap, in which case 

/duj ( J \ ^ 

— S inc (q,u) = a(-\ aq 2 , qa < 1, (22) 

Note that Eqs. (21) and (23) match each other (up to a 
numerical factor) at the transition line in Eq. (7). Com- 
paring Eqs. (18) and (23) with Eqs. (17) and (22), re- 
spectively, we make two interesting observations. First, 
at J/U -C (N/Nq) 2 the responses of the compressible 
and incompressible parts of the cloud occur at well- 
separated intervals of oj. Second, the total weight of 
the response of the Mott phase is smaller by the fac- 
tor Tr(N/N ) 2 (J/U)(aq) than the weight provided by the 
compressible phases, although the latter occupy only a 
small fraction of the atomic cloud. 

The decrease of the integral intensity in the Mott phase 
does match the experimental observations, and does not 
contradict sum rules, see Sec. VB. Indeed the major 
part of the spectral weight in the presence of the lattice 
is shifted to higher frequencies associated with the tran- 
sitions from the lowest to higher single-particle bands. 



III. DYNAMIC STRUCTURE FACTOR OF THE 
HOMOGENEOUS MOTT PHASE 



Hamiltonian into a main part (-ffo) and a perturbation 
(Hi) as follows 

#bh = H + H u (24) 

h o = !e(^<-'?) 2 ' (25) 

i 

H x = -jY J {b\ +1 h+b\b l+1 ). (26) 
i 

Here, rj controls the occupancy of a lattice site and is 
related to the chemical potential p in a linear fashion, 
j) = 1/2 + n/U. 

Upon raising the chemical potential p the lattice be- 
gins to be occupied at p = — 2J, and up to fj, = 
2 J [1 - 2J/U + O (( J/U) 2 )] , the state of the lattice is 
compressible, turning as next into a Mott-insulator state 
with site occupancy one. The Mott-insulator state turns 
to a compressible state at p, = U — 4J [1 + O (J/U)] and 
emerges again at p = U + 4J [1 + O (J/U)] with site oc- 
cupancy two. We shall be interested, in this section, in 
any Mott-insulator state, which amounts to choosing val- 
ues of 77 away from half-integers by a value ~ J/U. More 
precisely, to leading order in J/U, the state of the cloud 
is a Mott insulator if 

|/i — uU\ > 2(v + l)J, (27) 

for all v = 0, 1, 2 ... . The integer v for which l.h.s. in 
Eq. (27) takes the smallest value, gives the site occupa- 
tion number, 

p = v+- + -sgn(ii-vU). (28) 

The dynamic structure factor for the Bosc-Hubbard 
model is defined as follows 

S(q, uj) = -J2 dre^-wl ( Pl ( T ) Po ) , (29) 
a j J-00 

where pi = b\bi— (bjbi) is the deviation of the occupation 
number of site I from its average value and pi(r) denotes 
the Heisenberg representation of pi. 

Before going into the details of a rigorous derivation 
of the peak structure of S(q,uj) at u ~ U, we outline 
first a simpler derivation which gives the correct spectral 
weight of the peak. The derivation of the peak shape is 
deferred to Section IIIB. 



In this section we consider the Mott phase of the Bosc- 
Hubbard mode and evaluate the dynamic response to 
leading order of perturbation theory in J/U 1. The 
trap potential V(x) causes no inhomogeneity of the Mott 
phase and affects only the length over which the Mott 
phase extends. Considering here a macroscopically large 
Mott phase (extending over many lattice sites), we omit 
the trap potential V(x) from Eq. (1) and separate the 



A. Spectral weight at u> w U 

For a Mott-insulator state, the ground state in pertur- 
bation theory reads 

\*o) = (l-^)\*o) + 0(H 2 1 ), (30) 
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where |\&o) stands for the unperturbed ground state 



I*, 



n 



(31) 



with p being the Mott state occupation number. Note 
that the form of Eq. (30) is owing to the fact that Hi 
acting on |\l/o) creates only one kind of excitations. The 
excitation energy is approximately U and each excita- 
tion is represented by a particle-hole pair with hard-core 
interaction between the particle and the hole. 

Since here we are interested only in the spectral weight 
of S(q 7 uj), it suffices to write down a complete set of 
states spanning the space of all excitations with energy 
u> as U in the lowest order of perturbation theory. Or- 
thogonality with respect to Eq. (30) implies 



l*W a > = l*W»> + l*o) 



U 



0(Hf), (32) 



where I^Zjjj) stands for an unperturbed excited state 
with a particle at site l\ and a hole at site I2 =/= h, 

|tt llJa ) = 7 =L= blb la \*o), (h + h). (33) 
\fP\P + 1) 

Evaluating first the matrix elements 



Jy/p(p+l) 



U 



(34) 



we obtain for the dynamic structure factor in Eq. (29) 
the following result 



S{q,(j) 



I6irp(p + 1) 



— J sin 2 (qa/2)S(uj - U) . 



(35) 

Note that the peak weight increases with p like p 2 and is 
suppressed at small J by the factor (J/U) 2 . By gener- 
alizing the above derivation, it becomes clear that there 
exist also peaks in S(q, to) at frequencies equal to integer 
multiples of U, i.e. w ~ jU, however the weight of those 
peaks is suppressed as (J/U) 2 ^ . 

The (5-function in Eq. (31) comes from the fact that 
we neglected the dispersion relation of the excitations, 
assuming that they all have energy u> = U. As a mat- 
ter of fact, accounting for the dispersion relation of the 
excitations results only in a broadening of the ^-peak, 
without changing its spectral weight, see next Section. 



B. Semi-circular peak at u) ss U 

The states in Eq. (32) span the space of all excitations 
at energy ui rts U. However, they are not eigenstates of 
the Hamiltonian (24), not even at zeroth order of J/U. 
Below, we build up the eigenstates of the Hamiltonian 
(24) by a suitable linear combination of the states (32). 



In order to do so, we first consider the two-body problem 
that represents the excitation. We solve this problem by 
going to the center-of-mass and relative motion coordi- 
nates. Then, having proper states at the zeroth order 
of J/U, i.e. analogs of Eq. (33), we repeat the steps in 
Sec. Ill A to obtain the dynamic structure factor. 

To simplify our notations in this section, we set the 
lattice spacing to unity (a = 1) and agree to denote any 
coordinate dependence by a subscript and momentum 
dependence by a superscript. Our starting point are the 
states in Eq. (33) and the matrix elements of the Hamil- 
tonian (24) in the space spanned by these states. The 
interaction part Hq is already diagonal 



(% 1P2 \H |tf Wa 



3pihfip 2 l 2 



u - 



-E 
2 ^ 

1 



(p v) 2 



(36) 

and moreover, it is proportional to the unity matrix, 
and thus, is invariant under unitary transformations. We 
leave Hq aside for the time being and focus on the tun- 
neling part 



PlP2l^l|*/l(2> = ~J[(P+ 1)<W 2 (<Wl + l +<Wl~l) 



2+1 



,h-l 



)]• (37) 



We decouple the motion of the center of mass [28] of the 
excitation by the following transformation 



I*; 



1 



<l/2 



(38) 



where r = ±1, ±2, ±3, . . . , the sum over I runs over all in- 
tegers of the same parity as r (hence the notation I 6 /,.), 
and yt G [— tt,tt]. The states in Eq. (38) are normalized 
as usually, 



(39) 



In the new basis, the center of mass motion decouples 
and the matrix becomes diagonal in the momentum yc, 

(^Iffil*?) = WrvCxWx-V). (40) 
The resulting relative-motion problem reads 

TL r 'r{x) = —TteSr'^-i — T*5 r ' fr+ i, (41) 

and the effective tunneling 



where r',r = ±1,±2,±3, 
amplitude is given by 



T 



J 



(pe- 1 "' 2 + (p+l)e^/ 2 ) . 



(42) 



The momentum x enters in Eq. (41) as a parameter. 
Note that the range of values of r',r in Eq. (41) does 
not include zero. We can, therefore, split H. r ' r into two 
semi-infinite chains, H r > r = Ti^, r + H~, r , where H^, r is 
given by Eq. (41) with solely positive/negative values of 
r' , r. Each of the semi-infinite chains can be solved easily 
by fictitiously including zero into the domain of values of 
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r', r, but requiring that the wave function vanishes at site 
r = 0. The solution then reads 



sin(Ar) \ &± r ) , 



(43) 



were k € (0, 7r) and the sign index ± stands for the 
right/left chain. The phase (j)^ originates from T x = 
\TJe^" with 



\T„\ = J^l + Ap{p+l)co&{x/2). 
The states in Eq. (43) are normalized as follows 

= S s , 3 5(k' - k)6(K' - x). 



(44) 



(45) 



Taking into account the part Hq, we obtain that the 
particle-hole excitation has the dispersion relation 



e±(k,x) = U -2\T Jt \cos{k), 



(46) 



with \Ty\ given in Eq. (44). 

The next step would be to correct the states j^* 3 *) in 
perturbation theory to the first order of J/U, and thus, 
obtain the analog of the states in Eq. (32). This can be 
done straightforwardly, and for reference, the result reads 



2JVMP+1) 



u 



sin(fc)5(x)|*o)- (47) 



It is, however, not necessary to repeat the same steps 
as in Sec. Ill A — remarkably, all those steps can be com- 
prised in a single formula and one can use the bare states 
1^^*) to evaluate matrix elements. We use the following 
expression, which turns out to be a considerable shortcut 
in the calculation, 



(*o|^|*^) = - 



(48) 



where J\ is the current operator Ji = i[H, b\b{\ which 
reads explicitly 

Jl = iJ (bjb l+1 - b\ +1 b^ + iJ (bjbt-t - ftf^ft,) . (49) 

Equation (48) holds for the considered model to first or- 
der in J/U. 

Using Eq. (48) or repeating the steps in Sec. Ill A, we 
finally obtain 



(^olbjhl^ 



± 2Jy/p(p+l) cTU>x+ix i 
iirU 

x sin(fc) sin(x/2). 



(50) 



With the help of the matrix elements in Eq. (50) and the 
dispersion relation in Eq. (46), the dynamic structure 
factor in Eq. (29) evaluates to 



S(q,w) 



2p(p+l) 
a 

1 



4,7 

Jj 



sin 2 (qa/2) 



W p {q) 



1- 



U 



W P (q) 



(51) 



Here, we have reintroduced the lattice spacing a and de- 
noted the width of the peak by 2W p (q), where 



W p (q) = 2J v /l + 4p(p+l)cos 2 (<?a/2). (52) 

This completes the derivation of Eq. (17), which is read- 
ily obtained from Eq. (51) by specifying the value p — 1 
for the Mott lobe occupancy. Note also that the spec- 
tral weight of S(q,uj) in Eq. (51) is equal to the spectral 
weight of S(q,u>) in Eq. (35), as expected from our dis- 
cussion in Sec. Ill A. 

Finally, we note that the same problem was studied 
in Ref. 29 in the limit p 3> 1 using a combination of 
perturbation theory and periodic boundary condition on 
the wave function. The result of Ref. 29 for S(q,co) rep- 
resents a series of closely spaced (5-spikes located in the 
neighborhood of oj = U. These spikes become increas- 
ingly denser and can be replaced by a continuous curve 
in the limit of large period of the boundary condition. In 
order to compare our result against the result of Ref. 29, 
we replace the sum in Eq. (9) of Ref. 29 by an integral 
(an extra factor of 1 /2 is needed to account for the sparse 
summation in the sum). After integration, we recover a 
semicircular peak in S(q, w) centered at w = U. The 
width of the peak agrees with our result for width in the 
limit pcos(qa/2) 3> 1; the prcfactor of S(q,u>) agrees as 
well, provided differences in the definition of S(q,u>) are 
taken into account. 



IV. AVERAGED DYNAMIC STRUCTURE 
FACTOR 

In this section, we study the effect of inhomogeneity of 
the atomic cloud on the dynamic structure factor S(q, w). 
The inhomogeneity is caused by the trap potential V{x), 
see Eq. (2). The characteristic length over which the 
atomic density n(x) varies along the trap is usually much 
larger than the quantum length scale A of the longitudi- 
nal trap confinement. For such "soft" confinements, we 
employ the density averaging in Eq. (12). We first find 
n(x) in the regimes of interest, using the Thomas- Fermi 
approximation, and then evaluate the dynamic structure 
factor averaged over the trap inhomogeneity. 



A. Atomic density profile in a ID trap 

Here, we consider a ID Bose gas subjected to the po- 
tentials of a harmonic trap and an optical lattice. We 
assume that the lattice potential is strong enough to jus- 
tify our use of the Bose-Hubbard model in Eq. (1). As 
a limiting case of the Bose-Hubbard model (continuum 
limit), we recover the model of ID Bose gas subjected to 
only a harmonic trap potential, sec end of this subsection. 

For a soft trap confinement, the density profile n(x) 
is commonly found using the Thomas-Fermi approxima- 
tion. The trap potential V(x) is treated classically and 
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taken into account through the local electro-chemical bal- 
ance 



d£_ 

dn 



+ V(x) = li, 



(53) 



where £ is the ground state energy (per unit volume) of 
a homogeneous system of density n. In calculating £, 
energy should be measured with respect to the chemical 
potential. The profile n(x) is found by solving Eq. (53) 
for n. This approximation yields n(x) accurately for the 
most part of the cloud, with exceptions being the phase 
boundaries (points at which na = 0, 1, 2 ... ), where the 
approximation works only qualitatively. For the purpose 
of calculating the dynamic response of the whole cloud 
Eq. (53) is sufficient, because only the gross features of 
n(x) matter. 

For weak interaction U / 1 Jn$a <C 1 (Bogoliubov limit), 
due to the bosonic nature of the problem, the resulting 
density profile is due to the balance between mean-field 
interaction energy and external harmonic trapping. It is 
convenient to measure energies from the bottom of the 
band, introducing p = u + 2J. The density profile is then 
given by 



n(x) 



ft - V(x) 
Ua 



n 



1 



L 



(54) 



where uq = fx/Ua is the atomic density at the trap cen- 
ter, and 2L is the cloud length, with L = a-u ft/t§. For 
the Bosc-Hubbard model (1), Eq. (54) is obtained from 
Eq. (53) in the limit J/U 1, using perturbation the- 
ory in the interaction. At smaller bandwidth, J/U < 1, 
Eq. (54) continues to describe the gross features of n(x), 
as long as U / Jn a <C 1, sec discussion in Sec. II A. Equa- 
tion of state (relation between p and N) is obtained from 
Eq. (54) by integration, i.e. imposing the normalization 
condition on the density profile, ./V = J_ L n(x)dx, yield- 
ing 



N/No = - (p/U) 



3/2 



(55) 



where Nq is defined in Eq. (3). Upon substitution of 
p = noaU in Eq. (55), we recover Eq. (4). 

In the opposite limiting case V /Jn^a 3> 1, the 
Thomas-Fermi approximation is equivalent to the orig- 
inal Thomas- Fermi approximation for fermions, except 
for the fact that several layers of atoms can be added 
on top of each other at the price of raising the chemical 
potential by about U for each layer. For example, the 
second layer of atoms begins to form in the center of the 
trap when the following condition is satisfied (retaining 
leading order in J/U <C 1) 



U - 4 J = e (N/2)" 



(56) 



where 4 J represents half of the width of the second Mott 
band (see Fig. 2) and N/2 represents L/a at the leading 
order of J/U. In terms of the total number of atoms N 



and in the limit J — > 0, the second layer starts forming 
with increasing N at N = 2Nq. In this section, we will 
consider not more than a single layer of atoms in the 
trap, which can be rephrased, with the help of Eq. (56), 
as follows 



N/N < 2y/l -4J/U w 2 (1 - 2 J/U) 



(57) 



Precisely this linear dependence of the maximal N/N 
on J/U is displayed in Fig. 1 at small J/U (see line ©). 
Note that Eq. (57) is accurate only to first order of J/U, 
because we used leading order expressions for the Mott 
bandwidth and atomic cloud length in Eq. (56). 

At J/U <C 1 (Tonks-Girardeau limit), the solution 
can be found by mapping the bosons onto noninteract- 
ing fermions, where the repulsive on-site interaction can 
be accounted for by the Pauli exclusion principle. The 
atomic density n(x) in the compressible region is then 
found from 



2 J [1 - cos (iran(x))] + V(x) = p, 



(58) 



where we have used the tight-binding expression for 
the chemical potential of a homogeneous system of 
non-interacting fermions. In the incompressible region, 
Eq. (58) does not hold and n(x) is trivially given by 
n(x) — l/a. From Eq. (58), we find the well-known ex- 
pression for the atomic density 



2 . L-L-V(X) 

n(x) = — arcsm \ . 

ira V 4J 



(59) 



As before, the equation of state is found by integrating 
n(x) over x. At < p < 4J, the atomic cloud consists 
of a single compressible phase (state C in Fig. 1). From 
Eq. (59), we obtain 



N = 



4 \ U -p 



E 



p,-4J 



K 



/i-4J 



(60) 



where K(m) and E{m) are the complete elliptic integrals 
of the first and second kind, respectively. At 4 J < p, < 
U — 4J, a Mott phase is present in the middle of the 
cloud, see Fig. 2. Adding the contributions of the Mott 
phase and two compressible caps to the total number of 
atoms, we obtain 



4 fn /4J\ 
N = -J—E — 
7r V eo \ ft J 



ft > 4J. 



(61) 



At the lower limit of applicability of Eq. (61), fi = 4J, 
a Mott phase nucleates in the center of the trap; with 
£7(1) = 1 in Eq. (61), one obtains Eq. (7) for N. 

In the low-filling regime p, 4J we recover the con- 
tinuum limit. An arbitrary interaction strength can be 
accounted for within the Lieb-Liniger model [3] for bosons 
with effective mass defined by optical lattice. In this 
limit, all atoms reside at the bottom of the band, noC 1, 
and one can approximate the tight-binding dispersion re- 
lation as follows 



e{k) = 2,7(1 -cos (to)) 



k 2 
2m* 



(62) 
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where m* = 1/2 J a 2 is the effective mass in the contin- 
uum model. In the weak interaction limit, n(x) is given 
by Eq. (54). At strong interactions, Uj ' Jn$a ^> 1, one 
obtains from the fermionic representation of the model 
that 



n(x) = n \/l - 



(63) 



where tlq = (1 / wa) y/ ft / J and L = ay ft/eg. In 
both cases, n(x) is quadratic in x around the trap 
center — a feature that also persists through the crossover 
regime. [30] The equation of state is obtained by integrat- 
ing n{x). At U / Jn^a 3> 1, integration of Eq. (63) yields 



N 



/' 



(64) 



Upon substitution ft = J (irriQa) into Eq. (64), we re- 
cover Eq. (5). We note that Eqs. (63) and (64) are ap- 
plicable at arbitrary ratios J/U, as long as U/Jn a ^> 1, 
which means ft/ J -C min {l, (J7/J) 2 }. 

Having discussed the continuum limit of lattice bosons, 
it is easy to adapt Eqs. (54), (55), (63), and (64) to the 
case of confined bosons without an optical lattice. It is 
achieved by the replacement: m* — > m and Ua — > g. 



B. Dynamic structure factor in the continuum 
limit, Bogoliubov regime 

For sufficiently small values of U, the interaction strength 
is small, U/Jn^a <C 1, and the chemical potential lies 
within the lower parabolic part of the spectrum in the 
Bose-Hubbard model, p -C 4J. Considering also small 
transition frequencies to -C 4J, we may take the contin- 
uum limit, a ~ \/\fj ~ 1/U ~ y/eo — > 0, in the Bosc- 
Hubbard model and match it with the model for the Bo- 
goliubov sound in ID, with the only difference being the 
presence of the trap potential V(x). The Hamiltonian 
(1) can thus be represented by 



H 



dx(jr (x) 



1 d 2 
2m* dx 2 



+ V(x) - p 



cj>(x) 



~ / dxcfJ {x)^ (x)<f>(x)<f>(x), 



(65) 



where 1/m* = 2 J a 2 and 4>(x) is a boson field opera- 
tor normalized as [(f>(x), c/)'(x')] = S(x — x'). The model 
in Eq. (65) coincides with the model we discussed in 
Sec. II B, up to obvious notational differences: m «-> m* 
and g <-> Ua. 

In the homogeneous case (V(x) — 0), the dynamic 
structure factor at small U is sharply peaked at the fre- 
quency of the Bogoliubov spectrum, see Eqs. (9) and (10). 
In order to understand the limit in which Eq. (10) is valid, 
one should imagine taking t/ — ^ while maintaining p 
(i.e. Una) constant. Then, Eq. (10) is exact. In practice, 
Eq. (10) is a fairly good approximation to the dynamic 



structure factor at finite U for small U/Jna <C 1, since a 
weak interaction only broadens the <5-function in Eq. (10) 
to a power-law singularity, see Sec. IV D. 

Next, we average Eq. (10) with the help of Eq. (12), 
using the Thomas-Fermi density n(x) in Eq. (54). Wc 
express the result through the parameters at the trap 
center and, to avoid doubling notations, we redefine now 
v as follows, v = y^n aU/m*. The atomic density at the 
trap center no is given in Eq. (4). We obtain Eq. (14), 
with m — > m* . 

As mentioned above, the Bogoliubov limit is the limit 
of small U. The requirement U/Jnoa 1 can be rewrit- 
ten using the crossover boundary in Eq. (6) as 



( j 



3/2 



(66) 



whereas the requirement fi <SC 4 J can be rewritten using 
Eq. (55) as 



7V< A 



3/2 



(67) 



The two requirements can be met only if J/U 1. Fur- 
thermore, in order to justify the use of Eq. (12), we have 
to require that the size of the cloud be much larger than 
the quantum length scale associated with the motion of 
a band atom in the trap potential. This imposes a condi- 
tion on N which becomes more rigid and replaces Eq. (66) 
at U < t 3 / 4 eo /4 , 



7V> 



1 



J 



3/4 



(68) 



We summarize here our result: the singular depen- 
dence of S(q,Lu) on u>, present at weak interactions in 
the ID model, is smeared out towards lower frequen- 
cies. Thus, the <5-function divergence present in Eq. (10) 
is weakened to a square root divergence as shown in 
Eq. (14). 



C. Dynamic structure factor in the continuum 
limit, Tonks-Girardeau regime 

Next we consider the limit of strong interaction, 
U/Jn^a ^> 1, and, for simplicity, restrict our consider- 
ation to small energies, p,oj <C 4 J. Following the same 
steps as in Sec. IV B, we arrive at the Hamiltonian (65) in 
the continuum limit. In the limit of strong interaction, 
the Hamiltonian (65) maps onto a free fermion model 
(Tonks-Girardeau gas) 



H = dx^(x) 



1 d 2 

2m* dx 2 



V{x) 



1>{x), (69) 



where ip(x) is a fermion field operator normalized as fol- 
lows {ip(x), tp^x')} = 6(x—x'). Note that the interaction 
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parameter U is not present in Eq. (69), because the ef- 
fect of the repulsive interaction is accounted for by the 
fermionic nature of ip{ x )- The correspondence between 
ip(x) and the bosonic field cj>(x) in Eq. (65) involves a 
non-local phase factor. However, the latter cancels out 
in quantities which do not involve permutation of parti- 
cles. In particular, the phase factor cancels out in the 
dynamic structure factor, since for the density operator 
we have <fi^ (x)4>(x) = ijv [x)ip(x). 

For a non-interacting model, the dynamic structure 
factor can be easily calculated for an arbitrary V(x). 
However, since we are interested here in the limit of weak 
trap confinement, wc employ nevertheless the density av- 
eraging in Eq. (12). In the homogeneous case (V(x) = 0), 
the dynamic structure factor is given by Eq. (8), with 
to* and Ttn 



reads 



(l/a)^pp. We 
oo and it 



to — > to" and nn = kp = y/2m*fi 
note that Eq. (8) is exact in the limit U / Jna 
is a fairly accurate approximation at any finite but large 
interaction strength U/ Jna 3> 1, provided the considered 
momenta are small q <C kp(U/Jna), see Sec. IVD. 

Averaging Eq. (8) with the help of Eq. (12) and using 
n(x) in Eq. (63), we find S(q,ui) and express it through 
the parameters at the trap center. To avoid doubling 
notations, we redefine now kp to refer to the trap center, 
kp = txjiq. The averaged dynamic structure factor is 
then given by Eq. (13) for q < 2imo- At larger momenta, 
q > 27T7io, a similar expression reads 



S(q,to) 



■y/(e+(q) -to) (uj-e-(q)) 



7m <r 

x6(e+(q)-ui)8(ui-e-(q)). 



(70) 



The Tonks-Girardeau limit is the limit of small den- 
sities. The requirement U / Jn^a 3> 1 can be rewritten 
using the crossover boundary in Eq. (6) as 



u 



3/2 



(71) 



whereas the requirement u <C 4J (or equivalcntly n$a <C 
1) can be rewritten using Eq. (64) as 



A< A 



1/2 



(72) 



This condition is compatible with the condition of hav- 
ing no Mott phase in the trap, see Eq. (7). Addition- 
ally, in order to justify our use of Eq. (12) in the Tonks- 
Girardeau limit we need to require that N 3> 1, which 
guarantees that the length of the atomic cloud is much 
larger than the trap quantum length. 

Next we analyze Eq. (13) in the limit of small mo- 
menta, q kp. In this limit the change due to the av- 
eraging is dramatic. The spectral weight of S{q,to), con- 
tained before the averaging within £-(<?) < to <£+(</), is 
nearly fully pushed out into the region < to < £-(?)■ 
We compare the spectral weights in the two regions in 
what follows. To leading order in q/kp <C 1, Eq. (13) 



-, . TO* / 2m*LO 

S(q,co) = —j2-- 

q V k F q 

for \u> — kpq/m*\ < q 2 /2m*, and 



ir> (73) 

kp 



S{q,u>) 



2m* z to 
klq 



k F q 



-1/2 



m to q 
2kp~ 



k F q 



k F q 



q 

2k F 



(74) 



for < to < kpq/m* — q 2 /2m* . The net spectral weight 
in Eq. (73) is given by 



e +dtO q [2q~ 

— S{q,to) = — J — . 
2ir dir V kp 



(75) 



To the same order in q/kp, the net spectral weight in 
Eq. (74) reads 



or \ q q [zq ^ 

— S( qi uj) = — -—J—. (76) 
2tt 2ir dir \ kF 



Comparing Eqs. (75) and (76) to each other, we find that 
the spectral weight remaining in the region e_ < lo < e + 
constitutes only a small part ~ y/ 1 qjkp of the total weight 
q/2n. Note that the total weight before averaging is also 
qj2ix. This behavior can be described on a simpler level 
as follows. 

Let us ignore the structure of the (rectangular) peak 
in Eq. (8) and replace it by a 5-function with the same 
weight, 



S(q,to) — qS(uj — kpq/m*). 



(77) 



This approximation can be rigorously justified in the 
limit q —> 0, because the peak width tends to zero faster 
(oc q 2 ) than the frequency at which the peak is centered 
(oc q). Carrying out the averaging for S(q,to) in Eq. (77), 
we obtain 



S(q,to) 



kp 



\J (k F q/m*) 2 - to 2 



(78) 



Equation (78) coincides with Eq. (74) to leading order 
in q/kp in the whole range of to, except for the neigh- 
borhood of to = kpq/m*. A more precise condition for 
the validity of Eq. (78) reads \to — k F q/m*\ > q 2 /m*. 
Equations (77) and (78) explain the gross features of the 
redistribution of the spectral weight present in Eqs. (73) 
and (74). Note that the square- root singularity present in 
Eq. (78) resembles the one obtained for the Bogoliubov 
limit in Eq. (14). Although the qualitative behaviors 
of Eqs. (14) and (78) arc similar, the Bogoliubov and 
Tonks-Girardeau limits can be distinguished from each 
other at small q by the w-dependence in the numerators 
of Eqs. (14) and (78). 
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We summarize here our results. The step-like be- 
havior of S(q,u>) at uj = e + (q) present in the homoge- 
neous system is smeared to a square-root dependence 
oc — uj, with £+(<?) evaluated at the density of the 
trap center. The step at lu = £-(q) is smeared in the 
same way (yjuu — e_) for q > 2kp. For q < 2kp, the 
smearing leads to non-zero spectral weight in the region 
< w < £-(<z); the behavior of S(q,ui) at lu = £-(q) 
remains non-analytic, with a diverging derivative on the 
lower frequency side. 



D. Dynamic structure factor in the continuum 
limit, at arbitrary interaction strength 

In the limits of weak and strong interaction, the dynamic 
structure factor S(q,u>) turns to zero if u> exceeds the 
Lieb-1 frequency e+(q;no), corresponding to the maxi- 
mal density no in the trap (see Sees. IV B and IV C). 
Contrary to this limiting behavior, at intermediate inter- 
action strength, S(q,cu) ^ at to > e + (q;n ). Such be- 
havior of S(q, lu) for trapped atoms reflects the evolution 
of S(q,Lu), with the interaction strength, in a homoge- 
neous system. 

In the Bogoliubov limit, a homogeneous system is char- 
acterized by a J-function response, see Eq. (10). Devia- 
tions from that weak-interaction limit leads to a replace- 
ment of the 5-function singularity by a power-law diver- 
gence [4, 5] with an exponent /xi < 1, see Eq. (11). The 
singularity remains nearly symmetric in lu — £+(q) and 
the exponent close to one, if the interaction parameter, 



mg/n, without optical lattice, 

m*Ua/n, with optical lattice, 



(79) 



is small, 7 <C 1. Upon increase of 7, an asymmetry 
develops in the peak structure (the higher frequency part, 
lu > £+((?), becomes suppressed), and jii is decreasing 
towards zero. At 7 ^> 1 (Tonks-Girardcau limit), the 
shape of S(q,Lu) approaches 9(e + — lu), see Eq. (8). 

In the following, we study how the singular behavior 
of S(q, lu) at Lieb-1 mode changes due to averaging over 
the density n(x) of the atomic cloud. Our basic assump- 
tion is that n(x) is a smooth function of x and can be 
approximated around its maximum at the trap center by 
a quadratic expansion 



n(x) &n + -n"(0)x 2 , 



(80) 



where uq = n (0). We note that this assumption about 
n(x) agrees with the numerical analysis of Ref. 30. 

The behavior of S(q,Lu) at Lieb-1 mode is summarized 
in Eq. (11), see also Refs. 4, 5. We find that the di- 
vergence in S(q,oj) at Lieb-1 mode is weakened due to 
averaging over the trap. We average S(q,ui) in Eq. (11) 
using n(x) in Eq. (80) and obtain S(q, lu) in Eq. (15) with 



the coefficients: 

A = 

B = 

C = 



27T/U1A 

sin(27r/zi 

7T/Lti^iA 
C0S(7T/Zl) ' 

v 1 K 



[C0S(7T^1) - V{\ , 



1 - 2/ii 

Here, the common factor A reads 



cm 
A = — 
qL 



Sen 



2Ae|n"(0)| 



1/2 



(81) 



(82) 



where c = y / 7r/F(l + ^ii)r(3/2 — ui) is a factor order 
unity (1 < c < 2.2) and Ae is a scale defined below 
Eq. (86). The nature of non-analytic behavior of S(q,Lu) 
is established by Eqs. (15), (81) and (82) for arbitrary 7. 
Next we provide a more detailed form of S(q,Lu), which 
can be matched to the results of Sections IV B and IV C. 



1. Limit of large q and arbitrary 7 
Here, we consider q qo, with 
9o = 



m v, 



7<1, 
7^1- 



(83) 



In this limit, the dynamic structure factor close to Lieb-1 
mode is given by [4, 5] (cf. Eq. (10)) 



S(q, 



Kg x ( ^-g+(g) 

e + {q) dl -" 1 I vq 



(84) 



where K = nn/m*v is the Luttinger liquid parameter, 
Hx = 1 - 1/2A", and S e (x) = (e/2) \x\ e ^ . Note that 
linie^o fie{%) = S(x). The Luttinger liquid parameter K 
is a function of 7 only and, for the Licb-Linigcr model, 
it has the asymptotes [3] 



K 



1 + 4/7, 



7 < 1, 
7 > 1. 



(85) 



The equality sign in Eq. (84) holds for the Bogoliubov 
limit 7 — > 0, whereas at finite 7 Eq. (84) gives S(q, lu) by 
order of magnitude in a finite range around Licb-1 mode, 
\lu — £+(<?) I < vq. Note that at 7 3> 1, parameter K — > 1 
and v —* nn/m*. 

Substituting Eq. (84) into Eq. (12) and using Eq. (80) 
we arrive at 



S(q, 



n vq 



Ae + L\j 2|n"(0)|Ae J Q ^ 



1 



W + £ ' 

(86) 

where Ae = n Q (de + /dn ) and W = (lu — e + )/vq. In 
obtaining Eq. (86), we linearized the n-dependence of 
the Lieb-1 mode, 

e+{q,n) « e+{q,n ) + (n-n ), 87 
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and retained n-depcndcncc only in the position of the 
singularity in Eq. (84). The upper limit of integration is 
taken to infinity, since keeping it finite goes beyond the 
accuracy of these approximations. Equation (86) gives 
the leading divergence (at /i > 1/2) in S(q, to) at W — ► 0. 
All parameters in Eq. (86) refer to the trap center. 
Performing the integration in Eq. (86), we obtain 



S(q,uo) 



Q I n vq \ 
4e+L\j 2|n"(0)|Ae | 

x r( Ml ) 

+6{co-£+) 



vq 



Ml -1/2 



1 



Co) 



Mi 



The last term in Eq. (88) is an additive constant and can 
be omitted at finite fix — 1/2 > 0. However, in the limit 
fii — > 1/2, which corresponds to 7 ^> 1, the last term 
plays an important role. In this limit, Eq. (88), as well 
as Eq. (15) with coefficients defined in Eq. (81) yield 



S(q,i 



In 



q 



imq/m* 



uj — e+ 



(89) 



In Fig. 5, we show the singular dependence of S(q,co) on 
co as given by Eq. (88). 

Next we compare Eq. (88) to Eq. (14) in the limit 
7 C 1. We use the Bogoliubov spectrum (9) for the the 
Lieb-1 mode and expand it for q 3> m*v, 



£+{q) ~ m*v 2 + 



2m* 



(90) 



Substituting here v = y naU /to* and differentiating 
£+(q) with respect to n, we obtain the characteristic 
scale Ae = n^oU = m*v 2 . Thus, one should expect that 
Eqs. (88) and (14) coincide at \co - e + \ <C m*v 2 . Taking 
the limit fi\ — ► 1 in Eq. (88) we obtain 



S{q,co) ~ nJ — 



1 



1/2 



z+-u>). (91) 



On the other hand, the same result is indeed obtained 
from Eq. (14), after expanding the r.h.s. in terms of 
co — e g , retaining the leading order term, and taking the 
limit q ^> m*v. 

We summarize here our result: the exponent fi\ of the 
power law divergence of S(q, co) at Lieb-1 mode, predicted 
in Refs. 4, 5, is reduced by 1/2 due to averaging over the 
trap. 



2. Limit of strong interaction 7> 1 and arbitrary q 

In this limit, the peak is asymmetric, [4, 5] 

in 



S{q,u) * 2L 



£+ — £- 



CO — £ + 

+(l + u)6(e+-u)] 



[(l-P)0(u-e+) 



(92) 



Here, < jii < 1/2 and < v < 1 are given by 



2 / 1 

/i! = — 1 arctanQ arctanQ, 

7T \ 7T 



1 



tan(^i 



(93) 
(94) 



with Q = q/m*Ua = irq/kp"/- Note that here we use a 
slightly different paramctrization of the peak asymmetry 
than in Refs. 4, 5. The Lieb modes at strong interac- 
tion assume the expression in Eq. (8). Equation (92) 
gives S(q,co) by order of magnitude in a frequency range 



around Lieb-1 mode, 



< 



The equal- 



ity sign in Eq. (92) holds for the Tonks-Girardeau limit 
(Q < 1), cf. Eq. (8). 

We proceed in the same fashion as in Sec. IV D 1, sub- 
stituting Eq. (92) into Eq. (12) and using Eq. (80), and 
arrive at 



S(q,co) 



m* n (e + - e_) 
2qL\] 2|n"(0)|Ae 

x [(1 - v)6{W + + (1 + v)6(-W - 0] 





1 


in 




w + z 





(95) 

where Ae = kpqjvn* and W — (lo — £ + )/(e + — e_). In 
Eq. (95), the parameters are taken at the trap center 
and the approximations made here are the same as in 
Sec. IV Dl. In contrast to the case of Sec. IV Dl, the 
integral here does not diverge at W = 0, if fj,i < 1/2. 
Nonetheless, Eq. (95) gives the leading order term of the 
non-analytic part of S(q, co) at Lieb-1 mode. The analytic 
part of S(q,co) is not known at this order of the asymp- 
totic expansion in W -C 1 ■ Performing the integration in 
Eq. (95), we obtain 



S(q,Lo) 



to* no(e+—e-) 



2qL y 2|n"(0)|Ae 



CO — £- 



tan {—) Q 



1-P 



(1 - D)9(co - £ + ) 

1 



(e+ - lo) 



Mi 



(96) 



We retained here an additive constant (last term) to en- 
sure proper behavior of Eq. (96) in the limit /ii — > 1/2, 
in which Eq. (96) takes the form of Eq. (89). In Fig. 7, 
we show the non-analytic behavior of S(q, co) of Eq. (96) 
at the Lieb-1 mode. 

At Q -C 1, Eq. (96) agrees with the Tonks-Girardeau 
limit, see Eq. (70). Indeed, substituting \i\ ps 2Q/tt and 
v w 1 in Eq. (96) , we obtain 



S(q,uo) 



V2m* 



CO — £j 



Ae 



1/2 



0(£+ - co). (97) 
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FIG. 7: Sketch of the dynamic response S(q, u)) as function of 
ui in the regime of "strong interaction and arbitrary q" (see 
Sec. IV D 2). The response has a non-analytic point (circle) 
at Lieb-1 mode £+(?) (solid line) and is strongly suppressed 
away from e+(q) at distances exceeding e+ — e_ (dotted line). 
The solid line is plotted using Eq. (96) with the choice Q — 1, 
corresponding to fii = 0.375 and v ~ 0.67; the dotted line 
illustrates the qualitative behavior of S(q,u>) away from the 
Lieb-1 mode. The scale of the ordinate axis is arbitrary. 



On the other hand, expanding S(q, uj) in Eq. (70) in terms 
of uj — £ + , we obtain the same result. 

To conclude, the singular behavior of S(q, w) at w = 
s+(q) is owing to the weak (quadratic) dependence of 
n(x) on x in the trap center. More generally, the expo- 
nent Hi of the singular dependence changes compared to 
the homogeneous case like \i\ — > fj,i — 1/d, where d is 
the power of the subleading order expansion of n(x) in 
terms of x around its maximum. Normally, n"(0) =/= 0, 
i.e. d = 2, which is the case assumed in Eq. (80). 



E. Dynamic structure factor of the compressible 
caps of the Mott phases 

In the Mott-insulator regime the atomic cloud con- 
sists of compressible and incompressible regions. Here, 
we evaluate the density-averaged dynamic response of a 
compressible region as well as its relative contribution to 
the total response mentioned in Sec. II C, see Eq. (18). 
The dynamic response of an incompressible region has 
been considered in Sec. III. 

We consider the simplest case illustrated in Fig. 2, 
where the center of the cloud is incompressible, with 
na = 1, and the ends are compressible, with n(x) given in 
Eq. (59). We assume J/U <C 1 and agree to measure the 
chemical potential /i from the middle of the first Mott 
band at x = 0. The length of one compressible region 
AL is then given by 



AL 

a 




O ((J/ M ) £ 



(98) 



We further assume that the trap potential is sufficiently 
weak, such that the compressible region is macroscopic, 
AL/a ^> 1, allowing us to use Eq. (12). 

In order to calculate the dynamic response, we map the 
Bose-Hubbard model in Eq. (1) onto a fermionic model 
in the strong interaction limit J/U <C 1 and at restricted 
site occupancy na < 1, obtaining 



H 



J £ (O + +£(*oi 2 -/i)/ l t /l, (99) 

I I 



where fi is a fcrmion annihilation operator, {/;,///} = 
5u>. As usually, // faithfully represents bi in the low en- 
ergy subspace (uj < U) up to a non-local phase factor. 
The phase factor, however, cancels out in the density op- 
erator, resulting in bjbi = ///;. 

In the homogeneous case (eo = 0), the dynamic struc- 
ture factor reads[25, 26] 



S(q,u) = 



a l T(q,uj) 



\J{Usin(qa/2)) 2 -oj 2 



(100) 



where 



F{q,u) = 6>(0 ++ )6»(— JIh ) + (9(fi_ + )0(-fi__), 

il SS ' = fi + 2sJcos (qa/2) cos^o + s'uj/2. (101) 

Here, s, s' = ±1 and 



arcsm 



4Jsin(qa/2) 



(102) 



It is convenient to relate averaging over the density 
profile in Eq. (59) to an averaging over chemical potential 
in expressions for the homogeneous case, 



S{q,u;ij) = 



2AL 



2.1 



S(q, uj; n')dfi' 
-2J y/eoifi- n') 



(103) 



Here, we assumed that 2 J < fj, < U — 4 J. 

Next, we average S(q,uj) in Eq. (100) using Eq. (103) 
and obtain that S(q,w) is given by r.h.s. of Eq. (100), 
but with T(q,uj) — > T(q,u>), where 



+ 2 J - Vm - 2 J 



(104) 



Thus, we have obtained the dynamic structure factor of 
a compressible region, normalized to its length. It is 
convenient next to expand Eq. (104) for fi 3> J, which 
corresponds to the compressible region occupying only a 
small portion of the cloud, AL <C L. To leading order, 
the dependence on [i drops out and we obtain 



S(q,u>) 



1 



2Ja 



\J (AJsm(qa/2)f - uj 2 



(105) 



with the support < uj < 4J| sin(ga/2)|. We note that 
this result could also be obtained by making a linear ap- 
proximation for V(x) over the size of the compressible re- 
gion, which is compatible with our assumption AL <C L. 
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Finally, we are in the position to obtain Eq. (18) for 
the contribution of one compressible region relative to 
the total length of the cloud. We need only to multiply 
Eq. (105) by AL/2L, where 



ifi + 2j nr 

a\l w a. —. 



(106) 



Expressing p via N with the help of Eq. (61), p = e iV 2 /4 
at p 3> J, wc obtain that AL/2L = AJ/eoN 2 and multi- 
plying it in Eq. (105) we arrive at Eq. (18). 

We remark that Eq. (104) holds also in the fully com- 
pressible regime — 2 J < p < 2 J, provided the square 
roots are assumed to vanish for negative arguments, i.e. 
redefining ^fx := 9(x)yfx. 



V. DISCUSSION 

A. Ways to measure dynamic structure factor 

Bragg spectroscopy[9-16] allows one to excite the 
atomic cloud with a scalar potential of the form 



ip(x, t) = Re ipoe 



i(qx-tvt) 



(107) 



with lu and q tunable independently of each other. The 
control of amplitude tpo, ideally, is independent of a; and 
q. The effective interaction of the atoms with the Bragg 
probe reads 



H int = / dxp(x)(p(x,t). 



(108) 



The linear response of the atomic cloud to such a pertur- 
bation is characterized by the dynamic structure factor 
studied in this paper. 

One way to measure the dynamic structure factor is to 
perform calorimctric measurements on the system after 
exposing it to the driving potential ip(x, t) for some time. 
The energy gained by the system per unit time is [31] 



— = u \ip Q \ 2 S N (q,uj) 



(109) 



where S^(q,u>) ~ 2LS(q,uj) is the (unnormalized) dy- 
namic structure factor of the system. 

In the experiment of Ref. 15, Bragg scattering was used 
to excite the atomic cloud at a fixed momentum q and a 
variable frequency uj for a fixed duration of time. The re- 
sponse of the atomic cloud to the excitation was observed 
as a smearing of interference pattern of the atomic den- 
sity at the final stage of experiment, when the atomic 
cloud was left to expand freely without confinement. Es- 
tablishing a relation between the smearing of the inter- 
ference pattern of the matter waves and the amount of 
absorbed energy goes beyond the purpose of this paper. 
Here, we discuss qualitatively only the limiting case, in 
which the system reaches thermal equilibrium before the 



trap confinement is released. In this case, one expects 
that the broadening of the interference pattern is propor- 
tional to the temperature of the system, which has been 
"heated up" by the Bragg excitation. The frequency and 
momenta at which the system was excited enter in the 
result only through Eq. (109), and thus, one may take 
the quantity tuSN(q,uj) as a measure for the calorimctric 
response. 



B. /-sum rules in the presence of optical lattice 

The response u>SN(q,u) features several properties 
which wc would like to mention: 

1. Complete f-sum rule. Regardless of the interaction 
strength, trap potential, and optical lattice, the particle 
number conservation for atoms imposes the constraint 



°° du Nq 2 
—u>S N (q,w) = - — . 
2tt 2m 



(110) 



where m is the bare atomic mass. According to Eq. (110), 
the net spectral weight of the calorimctric response is 
conserved and can change only with particle number N 
and probe momentum q. 

2. Partial f-sum rule I. A similar sum rule to the one 
in Eq. (110) can be formulated in the continuum limit of 
the optical lattice problem. For the cases considered in 
Sec. IV B, IV C, and IV D, the partial /-sum rule reads 



f 4JdUJ Q ( \ 

/ — ub N {q,u)) 
Jo 27T 



Nq 2 
2m* 



(111) 



Equation (111) is valid for q <C ft /a and provided that 
p, U -C 4 J. This allows us to choose the upper limit 
of integration to be on the order of the band width 
4 J and to use the effective mass m*. The remaining 
spectral weight oc (1/m — 1/m*) is distributed between 
higher energy bands, at frequencies above the recoil en- 
ergy Eji = n 2 /2ma 2 . The interaction affects the shape 
of the integrand in Eq. (Ill), but does not affect the net 
weight as long as interband matrix elements of the inter- 
action can be neglected. One may check that Sn (q, ui) 
evaluated with the help of Eqs. (13), (14), and (70) sat- 
isfies the rule in Eq. (111). 

3. Partial f-sum rule II. Within the Bose-Hubbard 
model, following the standard derivation of /-sum 
rule, [31] one finds 

/•OO J 

/ — uS N {q,u) = -2sin 2 (ga/2) {H x ) (112) 
Jo 2tt 

where H\ is the tunneling Hamiltonian in Eq. (26) and 
average (...) is taken over the ground state. The upper 
limit of integration in Eq. (112) is formally infinity for 
the Bose-Hubbard model. For an optical lattice in the 
tight-binding limit J, U <C Er, the upper limit of inte- 
gration can be replaced by Er. Equation (112) can be 
considered as a generalization of Eq. (Ill), to which it 



17 



reduces in the limit q <C it /a and p,, U <C 4 J. Wc remark 
that the following inequality — (ifi) < 2 JA^ holds for ar- 
bitrary occupation na < 1 in the Bosc-Hubbard model. 
Therefore, the spectral weight is bounded from above as 
follows 



Em j 

—u;S N (q,uj) <4JNsm 2 (qa/2). (113) 

Z7T 



The latter equation shows that the dynamic response in 
the low frequency range ui G [0, Er] is suppressed at least 
as a J. A suppression that is only proportional to J 
(or 1/m*) is due to the optical lattice and is present 
also for non- interacting atoms. The interaction leads to 
an additional suppression when entering into the Mott 
regime. 

4- Suppression of dynamic response in the Mott regime. 
In the Mott phase the tunneling is quenched to first order 
of J. i.e. (Hi) = 0(J 2 ). At the second order of J, wc 
find (Hi) = 2NA£q, where A£q < is the correction 
to the ground state energy per atom due to tunneling. 
Thus, for the Mott phase we have according to Eq. (112) 



pu;S N (q,Lu) = A\A£ Q \Nsm 2 (qa/2). (114) 
Zn 



Evaluating A£q by perturbation theory, we find A£q = 
~2(p+l)J 2 /U, where we kept for generality an arbitrary 
Mott occupation number p. One may check that Sn (<?, u>) 
evaluated with the help of Eq. (51) satisfies the rule in 
Eq. (114). Since |A£ | <C J, a natural question arises: 
How does the suppression of dynamic response occur with 
entering into the Mott regime? To answer this question, 
it is instructive to consider a homogeneous compressible 
liquid in the limit J/U <C 1. Considering na < 1, wc 
obtain from Eq. (112) 



Er dco „ , „ 9 , , „ sm.( wna) 

— wS N {q,u) = 4JNsm 2 (qa/2) — i '-. (115) 

2tt ima 



Here, the factor sin(7ma)/7ma, which vanishes at na — > 
1, has at its origin the phase space restriction imposed 
by the interaction. This factor governs the transition 
between the upper bound in Eq. (113) and Eq. (114) 
with entering into the Mott regime. Finally, we remark 
that Eq. (115) is generalized to larger occupations p < 
na < p + 1 by replacing J — > (p + 1) J and na — > na — p, 
where p is a positive integer. 



C. Relation to the experiment of Ref. 15 

In this section, we turn to a discussion of the dy- 
namic response measurement reported in Rcf. 15. At 
first, we consider the case without optical lattice, namely 
the case of Fig. 1(c) of Ref. 15. The dynamic response 
is a single broad peak centered at Ld/2n ss 5 kHz with a 
FWHM of about 7 kHz, see Fig. 8. The response comes 
from an array of ID tubes loaded with Rb atoms, with 
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FIG. 8: Attempted comparison of a theory for a single tube at 
T = with experiment for the case without an optical lattice. 
Experiment: The data points (x) are taken from Fig. 1(c) 
of Ref. 15. Theory: Several limiting cases are shown. The 
solid line shows the dynamic response in the Tonks-Girardeau 
limit (7 — > 00) given by Eq. (13). The dashed line shows 
the response evaluated in the Bogoliubov limit, see Eq. (14), 
at the value of mv 2 ~ 3.2 kHz. The dotted line shows the 
response of non-interacting bosons (7 = 0) given by Eq. (A3) . 



N « 250 atoms in the central tube and with a longi- 
tudinal trap frequency of lvq/2t: « 43 Hz. The system 
is probed at q m 7.3 x 10 6 m _1 . In Fig. 8, wc compare 
the predictions of our theory against the experimental 
data. The solid curve shows the dynamic response in 
the Tonks-Girardeau limit (see Sec. IV C) and is plotted 
using Eq. (13) with wtiq = kp = \f2~mp, and p, — Nujq. 
The dashed curve in Fig. 8 shows the dynamic response 
in the Bogoliubov limit (see Sec. IV B) and is plotted 
using Eq. (14) with mv 2 w 3.2 kHz. We explain below 
how we determine the value of the interaction energy 
mv 2 in this limit. The dotted line in Fig. 8 shows the 
dynamic response in the limit of non-interacting bosons 
(see Appendix A) and is plotted using Eq. (A3). Unlike 
the previous two results, the result shown by the dotted 
line takes into account the "quantumness" of the trap 
potential, see Appendix A. Note that we used no fit pa- 
rameters in Fig. 8, except for rescaling all curves to have 
the same area to obey the /-sum rule. 

The value of interaction energy mv 2 in the Bogoliubov 
limit was determined as follows. Using expressions for 
weak interaction, we relate p, to N and g, obtaining 



/i = mv 2 



4V2 



mgtuoN 



2/3 



(116) 



The coupling constant g can be extracted from the ex- 
perimental data in the Mott regime. In particular, we 
consider Figs. l(f)-(h) of Ref. 15 and single out a series of 
peaks which do not shift towards lower energies with in- 
creasing the height of the optical lattice potential. These 
peaks are attributed to the dynamic response occurring 
at multiples of U in the Mott regime. In the case of 
Ref. 15, we were able to identify the peaks occurring at 
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FIG. 9: The same as in Fig. 8, but with averaging over a 2D 
array of ID bosonic systems with different particle number 
iV, see Appendix B. The solid line (Tonks-Girardeau limit) is 
plotted using Eq. (B7) and the dashed line (Bogoliubov limit) 
using Eq. (B5). The dotted line (non-interacting bosons) is 
insensitive to averaging over N at zero temperature. 



u)/2tt « 2 kHz and 4 kHz with the expected dynamic re- 
sponse at U and 2C/, respectively. Our identification of 
peaks agrees with the one reported in Ref. 15. Having 
the value of U/2tt = 2 kHz at a height of the lattice po- 
tential of approximately 13i?R (see Fig. 1(g) of Ref. 15), 
we evaluate g using tight-binding theory, 



-s-^Ua, 

7T 



(117) 



where s is the potential height (s = 13) measured in units 
of En = it 2 /2ma 2 . The lattice potential is assumed to 
be V(x) = sEr sin 2 (ttx/o) and the lattice constant in 
the experiment is a = 415 nm. With this information, 
we obtain the value of mv 2 f=a 3.2 kHz. We note that 
this value of interaction energy corresponds to a value of 
the dimensionless interaction strength of 7 ~ 0.32. This 
agrees by order of magnitude with the value of 7 ~ 0.6 
reported in Ref. 15. The discrepancy may stem from a 
different definition of 7, which involves averaging over an 
ensemble of ID systems. [32] Our definition of 7 refers to 
a single ID system and it is as follows 



7 



mg 
n ' 



where no is the density in the trap center, 
definition, one can estimate 7 as 

7 = min{7<,7>}, 

4V2 mg 2 
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where 7< (7>) is 7 obtained using the expression for no 
in the Bogoliubov (Tonks-Girardeau) limit. 



FIG. 10: Dynamic response of free bosons at a finite temper- 
ature, see Appendix C. Data points (x) are the same as in 
Figs. 8 and 9 and are taken from Ref. 15. Solid line shows 
the dynamic response of non-interacting bosons in a single 
ID trap at T = 4.2 kHz (least-squeares fit) and N = 250. 
Dashed line shows the response of an ensemble of such ID 
systems with N distributed according to Eq. (Bl) and with 
N max = 250. Temperature is kept at T w 4.2 kHz to illus- 
trate the broadening of the peak due to ensemble averaging. 
Dotted line (close to solid line) shows the ensemble-averaged 
response at T ~ 3.3kHz (least-squeares fit). All curves are 
rescaled to have equal areas in accordance with the /-sum 
rule. 



Next, we return to the comparison made in Fig. 8. 
We note that the measurement data in Fig. 8 display a 
much broader maximum than the one predicted by the 
weak interaction (dashed) or free-boson (dotted) theory. 
The maximum predicted by the strong-interaction the- 
ory (solid) is sufficiently broad, but it is offset by a large 
energy (w 5 kHz) away from the position of the measured 
data. In order to gain an insight into the origin of the 
discrepancy we verify two hypotheses: (i) averaging over 
the array of ID traps is responsible for the broadening 
of the peak and (ii) finite temperature is responsible for 
the broadening of the peak. We carry out the averaging 
over the 2D array of ID systems in Appendix B. We take 
into account the fact that different ID systems may have 
different particle number N, using the probability distri- 
bution P(N) in Eq. (Bl) taken form Ref. 2. The distri- 
bution P(N) depends only on the particle number N max 
in the central tube of the 2D array. We take N max = 250 
and plot our results of Eqs. (B7) and (B5) in Fig. 9. Aver- 
aging over the 2D array modifies the analytic properties 
of Sn{q_,<^>) at the Licb modes, see Appendix B. The 
peak shape has changed significantly for the Bogoliubov 
and Tonks-Girardeau limits as compared to the case of 
a single ID system. Nonetheless, the width of the peaks 
is unaffected by the averaging. This can be easily under- 
stood by recalling that the ID system with the largest 
particle number N has the largest dynamic response (see 
/-sum rules in Sec. VB). Thus, the dynamic response in 
a 2D array of tubes is dominated by the tubes close to 
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the central tube, which have about the same number of 
atoms N ~ iV max . We conclude that averaging over N 
in the 2D array cannot explain the discrepancy seen in 
Fig. 8. 

Secondly, we verify hypothesis (ii) and analyze the ef- 
fect of finite temperature in Appendix C. 

We find that already noninteracting bosons at high 
temperature would yield a spectral function with a width 
of the order to the one measured in the experiments, see 
Fig. 10. This fact is not surprising, since the momentum 
q at which the system is probed is sufficiently large (a) 
to cause a sizable broadening of the free-boson peak (see 
dotted line in Fig. 8) and (b) to increase the exponent 
Hi of the power-law singularity of S(q, lo) at Lieb-1 mode 
in the homogeneous Licb-Liniger model. Despite the fact 
that effect (b) refers to the homogeneous case, the phys- 
ical picture behind it applies to any system: interaction 
effects are less pronounced at large momenta of the probe. 
In particular, for the Lieb-Liniger model at 7 ^ the ex- 
ponent /^i increases monotonically from fii = at q = 
to /Xi = 1 — 1/2K at q 3> qo, where go ~ max (mi), mj) 
is the momentum scale defined in Eq. (83). At q qo, 
a quantum (q, lo) of the probe excites single bosons pro- 
jectively rather than exciting many bosons collectively 
via processes involving interaction. For reference, we 
note [33] that at 7 = 0.6 the exponent Hi w 0.2 at 
q = O.I7TO0 and Hi f=a 0.8 at q = irn.Q. Similarly, effect 
(a) becomes pronounced with increasing q; the width of 
the peak scales as Slo ~ q/mX, see Appendix A, making 
any interaction-induced broadening difficult to observe 
at large momenta (q/mX 3> mv 2 ). A finite tempera- 
ture complicates the distinction between the interaction 
and quantum finite-size effects. In Fig. 10, the solid line 
shows the fit to the experimental data of the free-boson 
theory for a single ID trap, see Appendix C; the best 
fit is obtained at T = 4.2 kHz. A characteristic sharp 
peak of width Slo ~ q/mX with long tails extending over 
Slo ~ q^jT /m is indicative of the regime of high temper- 
atures, luq <C T -C loqN/IsiN. In the very-high temper- 
ature regime, T S> oo N/\nN, the sharp peak is absent 
and the dynamic response is a single broad peak. In 
Fig. 10, the temperature is intermediate between "high" 
and "very high", T ~ looN/ \nN, which explains the pe- 
culiar shape of the peak. In contrast to the case of T = 
(cf. dotted line in Fig. 9), ensemble averaging affects the 
free-boson dynamic response at T 7^ 0. Using the proba- 
bility distribution in Eq. (Bl), we carry out the averaging 
over N and show the result in Fig. 10 (dashed line); note 
that T for the dashed line was kept at the same value as 
for the solid line. Not unlike the results of Fig. 9, ensem- 
ble averaging does not lead to substantial broadening of 
the peak (width increses by factor order unity). Taking 
into account the ensemble averaging in our fit, we obtain 
a smaller value of temperature, T S3 3.3 kHz, see dotted 
line in Fig. 10. Since we neglected interactions, this value 
of temperature is an upper bound on the temperature in 
the experiment. 

In the following, we briefly discuss the case with an 



optical lattice. At sufficiently strong binding of the opti- 
cal lattice, the atomic cloud enters into the Mott regime, 
see Figs. l(f)-(h) of Ref. 15. As mentioned above, we can 
identify here several peaks which do not shift to lower 
frequencies with increasing the optical lattice amplitude. 
These peaks, in particular the two peaks occurring at fre- 
quencies S3 2 kHz and S3 4 kHz, can be attributed to exci- 
tations of the Mott phase occurring at frequencies U and 
2U, respectively At even lower frequencies, the peaks 
are most probably associated with the response of com- 
pressible domains, see Sec. IV E, at u < 4(p + 1) J, with 
p being the underlying occupation number, see Eq. (28). 
At low temperatures (T <C U), the main contribution to 
the peak at u> sa U comes from the incompressible phases. 
According to Eqs. (17) and (18), the contribution of the 
compressible phases at lo ~ J and incompressible ones at 
lo ss U into the /-sum rule are of the same order [34] at 
q ~ ir/a and N ~ -/Vo- This is in a qualitative agreement 
with the experiment, in spite of the fact that probably 
T ~ U in the experiment. The spectral weight of the 
peak at a; ss 2U is by ~ (J/U) 2 smaller than that of the 
lo S3 U peak, see Eq. (19). The tendency of the decrease 
of lo S3 2/7 peak weight is seen in Figs. 1(f) and 1(g) of 
Ref. 15. 

The dynamic response measured in Ref. 15 features a 
strong suppression in the range of frequencies < lo < 
Eji as the system enters into the Mott regime. Such 
a suppression is to be expected from the /-sum rule in 
Eq. (112). 



Conclusion 

We considered the effect of a smooth trap potential 
on the momentum-resolved inelastic light scattering off 
a one-dimensional interacting atomic gas. Cases with 
and without optical lattice superimposed on the trap po- 
tential are analysed in detail. The singularities in the 
scattering cross-section associated with the excitation of 
collective modes persist, though in a modified form, even 
upon the averaging over the varying along the trap gas 
density. Rounding of the singularities due to the confine- 
ment occurs only at the smallest momentum transfers q 
of the order q ~ l/^mu)o, where m and loq are, respec- 
tively, the mass and oscillation frequency for an atom in a 
trap. Another source of smearing, apparently dominant 
in the current experiments [15], is a finite temperature of 
the gas. 
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APPENDIX A: 



QUANTUM FINITE-SIZE 
EFFECT 



In the bulk of the paper, we consider the effect of a 
varying atomic density n(x) on the dynamic response by 
employing Eq. (12). There exists a competing effect to 
the one accounted for by Eq. (12). The competing effect 
stems from the finite size of the atomic could and it is 
most vividly seen in the non- interacting limit. In that 
limit, the averaging in Eq. (12) has no effect, because 
the Lieb-1 mode is trivially equal to £+(<?) = q 2 /2m and 
thus is independent of n. In what follows, we consider 
the example of non-interacting bosons and determine the 
magnitude of the competing effect. 

For non-interacting bosons, the dynamic structure fac- 
tor can be calculated directly, 



S N (q,cu) = 2nN 



(<?A) 



2/' 



1=1 



2 l l\ 



■5(cj-two), (Al) 



where A = 1/ y/mwo is the quantum length scale of the 
trap and loq is the trap frequency related to eo in Eq. (2) 
by ujq = (l/a)y t /2eo/m. For simplicity, we consider the 
case without an optical lattice. It is convenient to re- 
place the sum over the discrete i5-peaks in Eq. (Al) by 
a continuous envelope function. Taking the continuum 
limit of £ in Eq. (Al), we obtain 



S N {q,uj) = N 



2tt 



wwq \2mu> 



j/u>o 



exp 



- q 2 /2m 

UJ 



(A2) 

Furthermore, since we are interested in small ujq, we may 
take the limit q 2 /mu>o ^> 1 in Eq. (A2) and obtain for 
the leading order asymptotic term 



S N (q,uj) = N\l exp 



2n 



(lo - q 2 /2m) 
Lu q 2 /m 



(A3) 



In the limit ojq — > 0, Eq. (A3) becomes Sjv(<7, w) = 
2ttN6(u! — q 2 /2m) as expected for the case without a 
trap. 

We summarize the previous paragraph as follows. The 
(5-peak in S{q,oS), occurring at frequency £+(<?) = q 2 /2m 
for free non-interacting bosons, is broadened by 5u> = 
q/mX, due to the "quantumness" of the trap potential. 
The width Sui can be rewritten in a more general way 



Soj 



dq 



Sq, 



(A4) 



which holds also in the presence of an optical lattice. 
Here, Sq is the momentum uncertainty in the ground 
state due to the zero-point motion in the trap. For the 



case of optical lattice, e+ (q) is the band dispersion rela- 
tion and 5q ~ 1/A*, where A* is the characteristic exten- 
sion of the ground state wave function. For luq <C 4J, 
the length A* is given by the expressions above, with 
in — ► to*. 

For interacting bosons, both the momentum uncer- 
tainty and the density variation affect S(q,oj). Equa- 
tion (Al) gives the limiting behavior of the dynamic 
structure factor for non-interacting bosons (or extremely 
weakly interacting, p, <C ojq) and Eq. (12) gives the lim- 
iting behavior for a nearly classical ("soft") trap poten- 
tial. The crossover between the two limiting cases has 
a governing parameter, which involves the interaction 
strength, curvature of the trap potential, and the mo- 
mentum q at which the system is probed. To determine 
this parameter, we add to Eq. (A4) a term representing 
the smearing due to the density averaging employed in 
Eq. (12), 



de+ de+ 
6uj = ——On + — — dq, 
an oq 



(A5) 



where Sn is on the order of no and represents the density 
variation in the trap. The governing parameter is ob- 
tained by taking the ratio of the two terms in Eq. (A5). 
As an example, let us consider the limit of weak inter- 
action and a sufficiently large q, such that we may use 
Eq. (90) for the Lieb-1 mode. We obtain from Eq. (A5) 
that 5u) = \± + q/m*X* . Thus, the density averaging pre- 
scription of Eq. (12) is valid for 



q <C /iTO*A* 



(A6) 



This condition is compatible with the requirement q S> 
1/A, see Sec. II B, if p, ^> lu^. Similar (but not identical) 
criteria as in Eq. (A6) can be worked out using Eq. (A5) 
for each limiting case in which £+(q) is known. 

A different (though related) finite-size effect is the dis- 
creteness of energy levels in the trap given by the trap fre- 
quency luq. This discreteness introduces a "coarse grain- 
ing" of the dynamic response along the frequency axis. 
In practice, ujq is negligibly small. In the main text of 
this paper, we assumed that both ujq and the tempera- 
ture are the smallest energy scales in the problem and 
neither of them is resolved on the scale of u>. 



APPENDIX B: AVERAGING OVER AN ARRAY 
OF ID SYSTEMS 

Here, we consider an ensemble of ID systems with dif- 
ferent number of atoms N. We average our results with 
the following probability distribution [2] 



P(N) 



1 



3 N 2/3 N 1 / 3 ' 



N < N n 



(Bl) 



where N mayi is the maximal number of atoms per ID 
system found in the ensemble. The probability distribu- 
tion in Eq. (Bl) describes a situation typically realized 
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in the experimental. 2. 15, 16] Namely, a large number 
(10 2 -10 3 ) of ID systems form out of a 3D Bose-Einstein 
condensate after an exponential ramp up of a 2D opti- 
cal lattice potential. The distribution of atoms over the 
array of ID systems is obtained as a snapshot of the orig- 
inal 3D density profile (assuming that the atoms had no 
time to migrate between the ID systems). For a 3D Bose- 
Einstein condensate in a spherically-symmetric harmonic 
trap, the number of atoms per ID system at position (i,j) 
in the 2D array is given by [2, 35] 



N, 



1,3 



2ttK 



5AT; 



max tf +j 2S 
3D 



3/2 



(B2) 



where -/V3D is the total number of atoms loaded into the 
3D trap. The distribution in Eq. (Bl) is readily obtained 
from Eq. (B2), after counting the number of ID systems 
with Nij = N for a large 2D array. 

We start with the Bogoliubov limit analyzed for a sin- 
gle ID system in Sec. IV B. We express the unnormalized 
dynamic structure factor Sj^(q, uj) = 2LS(q,uj) explicitly 
though the number of atoms N in the ID sys tem. Us- 
ing Eqs. (4), (14), and (55) as well as L = a^/pjeo and 
£o = muj^a 2 /2, we obtain 



<SW(q,^) 



2V2tt uj 2 - q 2 /2m 



^0.9 q \Je+{q) -uj 2 



(B3) 



for q 2 /2m < uj < e + (q), and Sn(q,^) — otherwise. 
The structure factor in Eq. (B3) depends on N only 
through the speed of sound v entering the expression for 
s+(q) in Eq. (9). The speed of sound v depends on N as 
follows 



V4\/2m 



1/3 



jyi/3 
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We average S]sr(q,uj) in Eq. (B3) over N with the prob- 
ability distribution in Eq. (Bl) and obtain 



<7 ( \ 25/2n 
b N {q,u) = -g-^ — 



q''v*LUog V 2m . 



(B5) 



for q 2 /2m < uj < £+(</), and §N(q,u>) = otherwise. 
Here, we replaced iV m ax — > N after performing the av- 
eraging. We conclude that averaging over an array of 
ID systems modifies the analytic behavior of SN(q,ui). 
The one-over-square- root singularity present in Sn (q, uj) 
before the averaging is replaced by a square-root non- 
analyticity after the averaging. The averaging does not 
change the width of the peak of <Sat(9, w) versus w, pro- 
vided we agree to compare the averaged response against 
the response of a ID system with N — N max . 

Next, we turn to the Tonks-Girardcau limit analyzed 
for a single ID system in Sec. IV C. Using Eqs. (5), (13), 



and (64), we obtain 

2 3/2 m l/2 



S N {q,uj) 



1/2 
quj < 



(lu - q 2 /2my 



N — 



2uiQq 2 /m 



<N- 



(LO + q 2 /2mY 
2ujoq 2 /m 



(B6) 



where the square roots arc redefined as ^fx := 6(x)yfx. 
Averaging the dynamic response in Eq. (B6) over N with 
the probability distribution in Eq. (Bl), we obtain 



S N (q,u>) 



2 3/2 m l/ 2iV l/2 



1/2 
quj ' 




q 2 /2m) 



2co Nq 2 /m 



-f 



(uj + q 2 j2m)' 
2u Nq 2 /m 



where the function f(x) is defined as follows 
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(B8) 



with 2Pi(a, (3, 7, z) being the hypergeometric function 
and the square root (re)defined as before. Not unlike the 
previous case, the averaging does not change the peak 
width in Sn(q,w) versus to. The non-analytic properties 
of S i at((7, uj) are altered by the averaging, making the non- 
analyticities at u = e±(q) less pronounced. We note that 
f(x) w (4/9)(l - xf/ 2 at x -> 1. 



APPENDIX C: DYNAMIC RESPONSE AT 
FINITE TEMPERATURES 

In the bulk of the paper, we were considering dynamic 
response at zero temperature. Many of the peculiar fea- 
tures of S(q, to) vanish gradually with increasing the tem- 
perature. The temperature can usually be regarded as 
the low energy cutoff limiting the observability of the 
singular and non-analytic behaviors of S(q,Lo). A gen- 
eralization of our theory to finite temperature is possi- 
ble, although it goes beyond the purpose of this paper. 
In this paper, our goal was to establish the regimes in 
which the power-law singularities of S(q,uj) derived for 
the Lieb-Liniger model can be observed in a setup with a 
trap. The finite temperature results we present below are 
needed for explaining the data of Ref. 15, see Sec. VC, 
and are thus elucidating the conditions realized in that 
experiment. 

We consider non-interacting bosons at finite tempera- 
ture T in a trap of frequency loq. The dynamic structure 
factor is obtained as follows, 

S N (q,uj) = 2tt ^2 / dxe iqx ip* e (x)(pi<(x) 

x(l + N t )N t ,8(u;-u>u>), (CI) 
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where fe(x) is the harmonic oscillator state, Ng is the av- 
erage number of particles in state £, and u>u' = (C~ ^')^>o- 
Within the grand-canonical ensemble, Ne obeys the Bosc- 
Einstein distribution, Nt = [exp ((£u>o — pt)/T) — 1] , 
with the chemical potential 11 determined from the equa- 
tion of state, 



5> 

1=0 



N, 



(C2) 



at a fixed value of T. At thermal equilibrium, it is suffi- 
cient to consider only absorption processes (u> > 0), since 
the probabilities of energy emission and energy absorp- 
tion are related to each other by the detailed balance 
equation 



S(q,-w) = e- u/T S(q,u). 



(C3) 



Assuming hereafter oj > and evaluating the matrix el- 
ements of e lqx in Eq. (CI) on the harmonic oscillator 
states, we arrive at (cf. Eq. (Al)) 
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(f +N e+U )N U ,(C5) 



where L^x) are Laguerre polynomials. Temperature af- 
fects the dynamic structure factor in Eq. (C4) through 
the shape of the function JV(£). At T = 0, this func- 
tion is constant, AT(l) = N, and we recover the result 
of Eq. (Al). In order to analyze the behavior of Af(£) 
at finite T, we rearrange the sum in Eq. (C5) as follows. 
We represent (1 + N( +l/ )N L , as the series (one may regard 
this series as a formal expansion in powers of variable £") 



{l + N e+ „)N„ = J2a p (e)C p , 



a P {£) 



(C6) 
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where £ = exp (—ojq/T). Then we carry out the summa- 
tion in Eq. (C5) for each term of the sum of Eq. (C6) 
separately, using the identity [36] 
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where Ig(x) is the modified Bessel function of the first 
kind and \z\ < 1. As a result we obtain, 
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Further, we analyze the function s(£); the dynamic struc- 
ture factor in Eq. (C4) is expressed through s(£) as 



S N (q,uj) 



two). 
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We take the continuum limit over £, assuming that ojq 
is small, and thus, only large £ = lu/luq are of interest. 
We are interested in frequencies oj ~ q 2 /m, which corre- 
sponds to £ ~ (qX) 2 . The Bessel function in Eq. (CIO) 
has thus to be evaluated at both large orders £ ^S> 1 
and proportionally large values of its argument. In order 
to do this, we use Dcbyc's asymptotic expansion of the 
Bessel function [3 7] 
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where Uk (t) are polynomials of order 3fc of t = 1/vl + z 2 . 
We refer the reader to Ref. 37 for the explicit form of 
Uk(t). Here, we consider only the leading-order term in 
Eq. (C12), which amounts to dropping the sum over k. 
Expanding, for consistency, also other terms in Eq. (CIO) 
for large £ ~ (qX) 2 3> 1, we obtain 
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Equation (C14) is our main result in this appendix. We 
analyze it in greater detail below. With the help of 
Eq. (C14), the dynamic structure factor can be expressed, 
in the continuum limit over to, as follows 



2?r 

S N (q,uj) = — s {uj/uj ) . 



(C15) 



The convergence of the sum in Eq. (Cll) is guar- 
anteed by the factor e PAI / T in the summand. The ra- 
tio ll/T < is found for given values of N and ujo/T 
from Eq. (C2). We consider here only the case N ^> 1. 
With increasing the temperature, fi/T monotonically in- 
creases by absolute value from fx/T = —1/N at T = to 
fx/T w - ]n(T/u N) at T > uj N. As a result, the higher 
the temperature the faster the sum in Eq. (Cll) con- 
verges. This is in contrast to the sum in Eq. (C5) which 
converges fastest at low temperatures. The presence of 
the factor tanh(pwo/2T) in the exponent in Eq. (C5) 
introduces an additional scale p ~ T/oj , which can be 
larger or smaller than the convergence scale p ~ T / (— li). 
As a result, we obtain several temperature regimes: 

(i) T <C oj 0i low temperatures, 

(ii) ujq -cTcT,, high temperatures, 

(iii) T* <C T, very high temperatures, 
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FIG. 11: Dynamic structure factor Sn(q,^) of non- 
interacting bosons at different temperatures. Dotted line: 
at low temperatures, T <C cjq, the factor Sn(q, w) dis- 
plays a Gaussian peak of width Slo ~ •JuJoq 2 /m centered 
at lu ~ q 2 /2m, see Eq. (A3). Solid line: at high tempera- 
tures, wo < T < iVoN/]nN, the factor SN{q,io) retains the 
low-temperature peak at a reduced intensity and develops a 
much broader (background) peak of width Sui ~ \/Tq 2 /m 
centered at the same position if q 2 /m 2> T. Dashed line: 
at very high temperatures, u^N/lnN <C T, only the broad 
peak is present. The width of this peak continues to increase 
with temperature, see text at Eq. (C21). The plot is made 
using Eqs. (C14) and (C15) for a system of N = 1000 parti- 
cles at temperatures T/loq — (dotted), 100 (solid), and 250 
(dashed); a value of {qX) 2 /2 = 100 was used. 



where T» is defined as the temperature at which /i = 
— luq. In the limit > 1, we obtain from Eq. (C2) 
that T» satisfies the transcendental equation N — 
(T^/ujq) In (T*/u)q). Thus, to leading logarithmic order 
(IniV > 1), we have T* =u N/]n(N). 

Let us first consider the zero-temperature limit, T <C 
u)q. Approximating e _Wo/T w and fi/T ks -1/JV in 
Eq. (C14), we arrive at 
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The sum in the parentheses in Eq. (C16) equals N at 
iV > 1. Using Eq. (C16) in Eq. (C15), we recover 
the zero-temperature result of Eq. (A3), up to lu/luq *-* 
(qX) 2 /2 in the prefactor, which is within the accuracy of 
the asymptotic expansion. 

In the high temperature regime, wo"CT<T,,wc split 
the sum in Eq. (C14) into two parts. In the first part, the 
summation index p runs from p = 1 to p ~ T/luq. Here, 
we approximate tanh (pwo/2T) ss poj /2T and replace 
the sum over p by an integral. In the second part, the 
summation can be carried out explicitly as in the case 
of Eq. (C16), but with the lower summation bound at 



~ T/u>o- As a result, we obtain 
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where c is a number order one, which can be fixed by re- 
quiring ds(£)/dc = 0. In practice, choosing a fixed value 
for c in the range between 1 and 5 yields accurate re- 
sults. The integral in Eq. (C18) can be expressed through 
special functions if desired. The qualitative behavior of 
the functions Si(£) and S2{£), in the considered tempera- 
ture regime, is as follows. The function Si(i) is, up to a 
prefactor, the zero-temperature result, see Eq. (A3). It 
contributes to Sjv(<7, w) with a relatively sharp peak at 
lu pa q 2 /2m. In contrast, the function Si(£) has a much 
broader peak centered at the same place as the peak of 
S2CO: provided q 2 /2m T. The dynamic structure fac- 
tor S^{q,Lu) is, therefore, a sharp peak with broad tails, 
see solid line in Fig. 11. As the temperature increases, the 
relative spectral weight shifts from the peak to the tails. 
The width of the sharp peak is given by Eq. (A4) with 
Sq ~ 1/A. The width of the broad peak is given by the 
same equation but with 5q ~ \J mT. The physical expla- 
nation for the coexistence of the sharp and broad peaks 
is the fact that in this temperature regime a macroscopic 
fraction of particles occupy the ground state of the har- 
monic potential. The sharp peak is due to the particles 
in the ground state and the broad peak is due to the par- 
ticles in the thermal tail of the Bose-Einstcin distribution 
function. 

In the "very high" temperature regime, T » T,, the 
contribution S2(£) can be neglected. Furthermore, the 
upper bound of integration in Eq. (C18) can be extended 
to infinity, since the factor e^ z ^° suppresses the inte- 
grand before the upper bound is reached. As a result, we 
obtain 
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If q 2 /2m ^> T, then the dynamic response has a peak at 
lu pa q 2 /2m with width 8u> ~ {q/rri)^J (— /i)m. A simple 
expression for \x can be given only if temperature is suffi- 
ciently high, T > lu N, sec text below Eq. (C15). Then, 
the width of the peak is 



(C21) 
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If q 2 /2m <C T, then Sn(q, w) is a monotonically decaying We show the behavior of Siy(q, ui) in different temper- 
function of u> on the scale of temperature. ature regimes in Fig. 11. 
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